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Preface 


This  work  uses  a  modified  version  of  the  Advanced  Beam  Experiment 
which  was  originally  developed  by  the  Vibrations  Branch  of  the  Flight 
Dynamics  Lab  at  Wright-Patterson  AFB.  The  modified  version  of  the 
Advanced  Beam  Experiment  was  performed  by  Capt  Jacques  in  the  AFIT  labs. 

The  goal  of  this  work  is  to  create  a  finite  element  model  of  the 
Advanced  Beam  Experiment  using  MSC/NASTRAN.  The  modal  includes  damping 
inherent  in  both  the  structure  and  the  actuators.  The  finite  element 
model  is  used  to  perform  a  sensitivity  analysis  on  the  selection  of  the 
weighing  parameters  Q  and  R  used  in  Linear  Quadratic  Regulator  theory. 
Sensitivity  analysis  is  also  performed  on  two  optimization  techniques. 
The  first  involves  minimizing  the  performance  index  J  and  structural 
mass,  while  the  second,  Onada's  formulation,  involves  minimizing  control 
mass  and  structural  mass.  Chapters  II  and  III  cover  the  structure 
modeled.  Chapter  IV  includes  a  brief  review  of  finite  elements  and 
describes  how  the  model  was  constructed.  Chapter  V  is  a  review  of  LQR 
theory  and  describes  the  state-space  formulation  of  the  problem.  For 
those  interested  in  only  the  results  of  the  sensitivity  analysis, 

Chapter  VI  should  provide  all  the  necessary  information. 

This  work  would  not  have  been  possible  without  the  help  of  many 
people.  The  previous  work  of  Capt  Jacques  and  Capt  Cristler  provided 
the  details  on  which  this  work  is  based.  Dr.  Liebst  spent  many  hours  of 
his  time  working  with  me  to  clarify  optimal  control  theory.  Capt  Gans, 
my  thesis  advisor,  provided  guidance  while  my  committee,  Lt  Col  Bagley 
and  Dr.  Spenny,  gave  valuable  comments  so  that  I  could  clarify  the  final 
work.  Mom  and  Dad  were  always  there  to  listen  whenever  I  needed 
encouragement  and  support.  And  finally,  I  would  like  to  thank  my  wife 
Joan.  She  always  looked  for  ways  to  help  and  tried  to  understand  what 
it  was  that  I  was  doing.  She  sat  through  my  defense  and,  when  it  was 
over,  celebrated  with  me. 

Steven  L.  Story 
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Abstract 


There  has  been  much  interest  during  the  past  decade  to  develop  and 
launch  large  space  structures.  The  high  cost  of  launching  material  into 
orbit  will  require  that  these  structures  be  assembled  in  space  using 
light  weight  elements  which  are  vulnerable  to  dynamic  excitations. 

Active  control  may  be  necessary  to  rapidly  attenuate  large  amplitude 
vibrations.  The  active  vibration  control  system  is  usually  designed 
after  the  structure  has  been  optimized.  The  integrated  design  of  the 
control  system  and  the  structure  may  provide  additional  weight  savings. 
This  thesis  presents  a  sensitivity  analysis  of  the  structure/control 
optimization  problem.  The  structure  used  is  an  aluminum  rectangular 
beam  with  proof  mass  actuators  mounted  on  the  free  end  and  a  structural 
dynamics  shaker  attached  at  the  midpoint.  A  finite  element  model  of  the 
structure  is  developed  using  MSC/NASTRAN.  Linear  Quadratic  Regulator 
theory  is  used  as  the  control  law  with  velocity  feedback.  Constant  and 
variable  values  of  Q  and  R  for  the  performance  index  are  used.  The 
variable  values  of  Q  and  R  are  selected  to  minimize  total  system  energy. 
Optimization  methods  examined  are;  first,  the  minimization  of  the 
performance  index  J  and  structural  weight;  second,  Onada's  formulation, 
which  minimizes  control  weight  and  structural  weight. 
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SENSITIVITY  ANALYSIS  FOR  THE  SELECTION  OF 


LINEAR  QUADRATIC  REGULATOR  WEIGHING  PARAMETERS  Q  AND  R 
FOR  ACTIVE  VIBRATION  SUPPRESSION 
OF  A  CANTILEVERED  BEAM 


I.  Introduction 


During  this  decade,  the  deployment  of  large  structures  in  space 
will  become  a  reality.  NASA  is  currently  developing  Space  Station 
Freedom  and  the  Department  of  Defense  is  interested  in  systems  which,  if 
deployed,  will  require  large  stable  platforms  in  space.  Because  of  the 
high  cost  of  launching  material  into  orbit  and  because  of  launch  size 
restrictions,  these  built-up  structures  will  be  assembled  in  space  usir  ’ 
light  weight  truss-like  elements.  The  truss-like  elements  will  provide 
a  framework  to  which  mission-related  functional  components  will  be 
attached.  Because  of  their  low  weight,  however,  these  structures  will 
be  highly  flexible  which  will  make  them  vulnerable  to  dynamic 
excitations  from  a  variety  of  sources. 

Control  Approaches  To  Vibration  Damping 

A  structure's  response  to  dynamic  excitations  is  governed  by  its 
mass,  damping,  and  stiffness  characteristics.  The  large  size  of  space 
structures  coupled  with  their  light  weight  results  in  many  low 
frequency,  lightly  damped,  closely  spaced  vibration  modes.  These 
vibration  modes  must  be  controlled  in  order  to  minimize  their  effect  on 
the  system.  Passive  damping,  obtained  by  methods  such  as  sophisticated 
shock  absorbers  and  visco-elastic  coatings,  is  a  partial  solution  to  the 


problem,  but  it  cannot  control  all  of  the  vibration  modes.  Active 
control  may  be  required  to  rapidly  attenuate  large  amplitude  vibration 
modes . 

The  early  history  of  active  control  research  for  Large  Space 
Structures  (LSS)  is  well  documented  in  Ref.  1.  The  authors  of  this 
paper  suggest  four  challenges  (1:515)  which  face  researchers  in  the 
years  ahead.  These  challenges  are: 

1.  Design  control  systems  sufficiently  robust  so  that 
errors  in  structural  modeling  can  be  accommodated. 

2.  Establish  reasonably  accurate  structural  models. 

3.  Develop  auxiliary  control  laws  which  adequately  reduce 
plant  excitations. 

4.  Establish  the  proper  choice  of  control  law,  sensors,  and 
actuators  to  maintain  the  shape  of  a  large  space  structure. 

To  address  these  challenges,  the  authors  divided  the  field  of  LSS 
research  into  structural  dynamics  and  control  theory. 

The  structural  engineer's  job  is  to  develop  a  simplified 
mathematical  abstraction  (model)  of  a  structure.  Since  all  of  the 
dynamic  characteristics  of  a  structure  cannot  be  modeled,  care  must  be 
taken  to  ensure  that  the  most  important  dynamic  characteristics  are 
accurately  represented  in  the  model.  A  fundamental  question,  however, 
must  be  asked  (2:4);  "What  is  the  purpose  of  the  model?"  Will  an 
"exact"  solution  be  required  from  a  continuum  model  or  will  a  lumped 
mass  model  be  adequate?  Once  the  model  is  developed,  verification  must 
be  done  to  ensure  that  it  accurately  represents  the  structure. 

Several  methods  of  model  verification  are  presented  in  Ref.  2. 

The  most  common  method  used  is  to  compare  the  model  dynamic 
characteristics  to  experimental  results  on  the  actual  structure.  The 
large  size  of  space  structures,  however,  makes  it  impossible  to  test  the 
actual  structure  before  it  is  deployed  in  space.  Scale  modeling  of  LSS 
will  be  necessary  if  actual  ground  tests  are  to  be  conducted.  Past 
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experience  (3:924)  with  finite  element  modeling  in  aircraft  design 
suggests  that  finite  element  modeling  would  be  a  useful  tool  to 
dynamically  model  LSS.  The  system  of  linear  differential  equations 
resulting  from  finite  element  analysis  is  readily  expressed  in  state- 
space  form  for  use  in  designing  the  active  control  system. 

The  active  control  system  for  vibration  damping  of  LSS  generally 
consists  of  sensors,  controllers,  and  actuators.  Prior  to  1980,  there 
were  very  few  LSS  active  control  experiments  (4:471).  The  first  major 
U.S.  Government  program,  Active  Control  of  Space  Structures  (ACOSS),  was 
started  in  1978  and  completed  in  1984.  In  late  1989,  the  number  of  LSS 
experiments  (4:472)  "seem  to  be  approaching  flood  level."  This  "flood 
level"  is  the  focus  of  Ref.  4  which  surveys  the  literature  published  on 
experimental  LSS  work  accomplished  during  1985-1989. 

According  to  Ref.  4,  the  experience  of  many  researchers  suggests 
that  the  actuators  play  the  dominant  role  in  determining  the  success  or 
failure  of  a  LSS  experiment.  The  type  of  actuators  used  in  LSS 
experiments  can  be  used  to  lump  the  experiments  into  one  of  two  distinct 
categories.  The  first  category  of  experiments  are  based  on  grounded 
actuators.  Since  it  is  impractical  to  ground  an  actual  LSS  actuator, 
the  second  category  consists  of  those  experiments  using  inertial 
actuators . 

Previous  Work  on  the  Advanced  Beam  Experiment  (ABE) 

The  original  concept  of  the  ABE  developed  by  the  Wright  Research 
and  Development  Center  Flight  Dynamics  Lab  was  to  use  four  inertial 
proof  mass  actuators  mounted  in  pairs  to  control  bending  in  two 
orthogonal  planes  and  torsion.  Cristler  (5)  developed  the  actuator 
controller  and  then  demonstrated  active  control  using  a  Linear  Quadratic 
Gaussian  (LGQ)  design  and  modal  suppression  techniques.  Breitfeller  (6) 
demonstrated  active  control  by  using  a  low  authority  controller  based  on 
root  perturbation  techniques  and  a  high  authority  controller  based  on  a 
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frequency-shaped  cost  function.  Both  of  these  experiments  were 
partially  successful. 

Jacques  (7)  set  out  to  resurrect  the  ABE  in  what  is  now  known  as 
the  modified  ABE.  In  the  modified  ABE,  only  two  of  the  original  four 
proof  mass  actuators  were  available.  This  limited  the  experiment  to 
controlling  only  XY-plane  bending  and  torsion.  To  provide  better 
control  over  second  and  third  mode  bending,  a  structural  dynamics  shaker 
was  added.  This  violated  the  original  intent  of  the  ABE  which  was  to 
use  only  inertial  actuators.  It  did,  however,  open  the  door  for  further 
research  on  the  ABE. 

Control /Structure  Optimization 

The  ABE  represents  the  traditional  approach  to  controls  research. 
The  control  engineer  is  given  a  structure  and  told  to  develop  an  optimal 
control  system  to  achieve  some  dynamic  response.  This  approach  may  not 
prove  adequate  with  LSS  because  of  the  high  cost  of  placing  mass  in 
orbit.  Since  the  size  of  LSS  requires  that  they  be  constructed  in 
orbit,  they  do  not  have  to  be  designed  to  withstand  large  launch  forces. 
An  immediate  cost  savings  can  be  achieved  by  designing  very  light 
structures.  Light  structures,  however,  are  susceptible  to  vibrations. 
Passive  damping  can  attenuate  some  of  the  vibrations,  but  active 
vibration  suppression  will  be  necessary.  To  obtain  maximum  system 
performance  for  minimum  cost,  an  integrated  approach  to  structural/ 
control  optimization  is  necessary.  Refs.  3  and  8  survey  some  of  the 
issues  concerning  the  integrated  optimization  of  structures  and 
controls.  Two  of  the  important  issues  in  integrated  design  according  to 
(8:55-56)  are,  "cross-sensitivity  information”  and  "the  choice  of 
objective  function." 
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Problem  Statement 


There  are  three  goals  for  the  work  presented  in  this  thesis.  The 
first  goal  is  to  develop  an  accurate  finite  element  model  of  the  ABE 
using  MSC/NASTRAN.  The  model  will  include  damping  for  both  the 
structure  and  the  actuators.  The  second  goal  is  to  examine  the  "cross¬ 
sensitivity"  of  the  structure  and  the  control  system  to  the  weighing 
parameters  Q  and  R  from  Linear  Quadratic  Regulator  (LQR)  theory. 
Finally,  two  different  objective  functions  used  in  the  integrated 
structure/control  optimization  problem  will  be  examined  by  using 
sensitivity  analysis. 
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II.  Structure 


The  Advanced  Beam  Experiment  (ABE)  configuration  was  originally 
developed  by  the  Wright  Research  and  Development  Center  Flight  Dynamics 
Lab  to  model  the  large  space  structure  characteristics  of  low  frequency, 
lightly  damped,  and  closely  spaced  vibration  modes.  The  structure  and 
its  properties  are  presented  in  this  chapter.  The  beam  theoretical 
equations  of  motion  for  xy-plane  bending  and  torsion  are  derived  and 
numerical  solutions  given. 

Structural  Properties 

The  beam  is  a  long,  solid  aluminum  beam  of  rectangular  cross- 
section.  It  is  suspended  in  a  vertical  position  with  a  circular  disk 
attached  to  the  free  end.  The  circular  disk  provides  a  surface  on  which 
to  mount  the  control  actuators  and  has  approximately  the  same  mass  as 
the  beam.  The  disk  provides  the  primary  component  of  rotary  inertia 
which  lowers  uhe  first  torsion  mode  frequency  so  that  it  is  within  the 
control  bandwidth  of  the  controllers.  The  beam  is  shown  in  Figure  2.1. 
The  beam  and  disk  physical  and  material  properties  are  given  in  Table 
2.1  and  Table  2.2  respectively. 

Equations  of  Motion 

The  equations  of  motion  for  the  beam  can  be  derived  using 
Hamilton's  principle.  Hamilton's  principle  can  be  stated  as  (9:199) 

t2  t2 

j6(T-Ve)dt  +  j5Wncc?t=0  (2.1) 

where 

T  =  total  kinetic  energy  of  the  system. 
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Figure  2.1  Advanced  Beam  Experiment  Configuration 
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Table  2.1  Beam  Properties  (7:26) 


Beam  Property  Description 

Value 

Units 

Beam  Length  (L) 

1.797 

m 

Y  Cross-Section  Width  (a) 

2.565  x  10'2 

m 

Z  Cross-Section  Width  (b) 

1.925  x  10'2 

m 

Cross-Section  Area  (A) 

4.539  x  ±u'4 

m2 

Young's  Modulus  (E) 

7.446  x  1010 

N/m2 

Shear  Modulus  (G) 

2.827  x  1010 

25 
' — 

3 

ro 

Beam  Density  (p) 

2.766  x  103 

kg/m3 

Beam  Mass  (m) 

2.455 

kg 

Y  Moment  of  Inertia  (Iy)* 

1.526  x  10‘8 

m4 

Z  Moment  of  Inertia  (Iz)* 

2.709  x  10‘8 

Torsional  Moment  of  Inertia  (J)* 

3.292  x  10'8 

m4 

* 

Polar  Moment  of  Inertia  (Ip) 

4.235  x  10'8 

m4 

*  See  Appendix  A  for  sample  calculations. 


Ve  =  potential  energy  of  the  system,  including  the  strain 
energy  and  the  potential  energy  of  the  conservative 
external  forces. 

5Wnc  =  virtual  work  done  by  nonconservative  forces,  including 
damping  forces  and  external  forces  not  accounted  for 
in  V. 


2-3 


Table  2.2  Disk  Properties  (7:26) 


Disk  Property  Description 

Value 

Units 

Disk  Diameter  (d) 

3.048  x  10’1 

m 

Disk  Thickness  (t) 

2.540  x  10'2 

m 

Disk  Mass  (md) 

4.986 

kg 

X  Mass  Moment  of  Inertia  (Idxx)* 

5.790  x  10'2 

kg-m2 

Y  Mass  Moment  of  Inertia  (Ijjyy)* 

2.895  x  10‘2 

kg-m2 

Z  Mass  Moment  of  Inertia  (Id  )* 

2.895  x  10'2 

kg-m2 

*  See  Appendix  A  for  sample  calculations. 


6 (  )  =  the  symbol  denoting  the  first  variation,  or  virtual 
change,  in  the  quantity  in  parentheses. 
t1f  t2  =  times  at  which  the  configuration  of  the  system  in  known. 

Equations  of  Motion  -  XY-Plane  Bending.  For  a  beam  in  bending, 
shear  deformation  will  occur  as  shown  in  Figure  2.2.  In  Figure  2.2, 
a(x,t)  is  the  rotation  of  the  cross  section  and  v(x,t)  is  the  total 
transverse  displacement  of  the  beam  neutral  axis  in  the  y  direction. 

The  shear  angle  P  is  defined  as 

<2.2, 

For  a  long,  thin  beam  undergoing  transverse  vibration,  however, 
the  equations  of  motion  can  be  approximated  by  using  the  Bernoulli  Euler 
assumptions  of  elementary  beam  theory.  These  assumptions  are  (9:193): 
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Figure  2.2  Beam  Shear  Deformation  (9:203^ 


-There  is  an  axis  of  the  beam  which  undergoes  no  extension  or 
contraction.  The  x-axis  is  located  along  this  neutral  axis. 

-Cross  sections  perpendicular  to  the  neutral  axis  in  the 
undeformed  beam  remain  plane  and  remain  perpendicular  to  the 
deformed  neutral  axis,  that  is,  transverse  shear  deformation 
is  neglected. 

-The  material  is  linearly  elastic  and  the  beam  is  homogeneous 
at  any  cross  section. 

-o  and  o  are  negligible  compared  to  o  . 

y  z  x 

-The  xy-plane  is  the  principal  plane. 
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In  addition  to  these  assumptions,  rotary  inertia  and  gravity  will  be 
neglected  in  the  theoretical  calculations.  Equation  2.2  becomes 


dv(x,t) 
dx 


a(x,t) 


(2.3) 


For  displacement  in  the  y-direction,  the  beam  internal  strain 
energy  is 


Ve=*fE(x)I,(x)(  ill )  2dx 


<Jx" 


(2.4) 


and  the  beam  kinetic  energy  due  to  translation  is 


T=^jp(x)A(x)(lX^lll)2dx 


FT 


(2.5) 


where 

-E(x)  is  Young's  Modulus. 

-Iz(x)  is  the  moment  of  inertia  about  the  z-axis. 
-p(x)  is  the  density. 

-A(x)  is  the  cross-sectional  area. 


For  transverse  loading,  the  virtual  work  is 


dWn 


|p(x, t)dv(x, t)dx 


(2.6) 


Substituting  equations  2.4,  2.5,  and  2.6  into  Hamilton's  equation, 
integrating  by  parts,  and  applying  geometric  boundary  conditions  gives 
the  equation  of  motion 


-ii[g(x)Xz(x)  d2v(x,t)  j+p(x)A(x) 

ax2  dx2 


a2v(x, t ) 


dtc 


=p(x, t) 


(2.7) 
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and  the  natural  boundary  conditions 


E(x)Iz(x)  |x=L=wz(t) 

dxd 

E(x)Iz(x)  a2lli£^|x<=Vy(t) 
dx3 

For  free  vibration,  equation  2.7  reduces  to 

-ii[F(x)Jz(x)  ]+p(x)A(x)i^-(iltii=0 

dx2  dx2  dt 2 

Assume  a  solution  of  the  form 


v(  x,  t)  =V  (x)  cos  ((i)t -a) 

Substituting  equation  2.11  into  equation  2.10 

2  2 

-?—(E(x)  Iz(x)  iLZip.]- p(x)A(x)V(x)a£=0 
dx2  dx 2 

For  a  uniform  beam 


dVoo  _pW  0 

dx4  ^ 


Let 


pAo^ 


then  the  eigenvalue  problem 


pV(x)=o 

dx4  v 


has  a  general  solution  of  the  form 


V(x)  --C1sin(  PyX)  +C2cos{  PyX)  +C3sinh(  PyX)  +C4cosh(  PyX) 


(2.8) 


(2.9) 


(2.10) 


(2.11) 


(2.12) 


(2.13) 


(2.14) 


(2.15) 


(2.16) 
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Equations  of  Motion  -  Torsion.  Hamilton's  principle  can  also  be 


used  to  derive  the  equations  of  motion  for  torsion.  The  potential 
energy  for  torsion  is  given  by  (9:200-202) 


L 


30 (x, t ) 
33? 


]2dx 


(2.17) 


and  the  kinetic  energy  is 

L 

r  =  ^[p(x)Jp(x)(il^_^i]2dx  (2.18) 


where 

-G(x)  is  the  shear  modulus. 

-J(x)  is  a  geometric  property  of  the  cross  section. 
-Ip(x)  is  the  polar  moment  of  inertia. 

-0(x,t)  is  the  rotation  at  x. 


The  virtual  work  of  the  external  forces  is 


3wnc=  |T  t  )30  (x,  t )  cfx+rL  (L,  t )  30  (L,  t ) 


(2.19) 


where 

-r(x,t)  is  a  distributed  moment. 
-TL(t)  is  a  concentrated  end  moment. 


Substituting  equations  2.17,  2.18.  and  2.19  into  Hamilton's  equation, 
integrating  by  parts,  and  applying  the  geometric  boundary  conditions, 
the  equation  of  motion  is  found  to  be 


p(x)Ip(x) 


d2d (x, t )  3 


3t‘ 


33? 


(G(x)J(x) 


30(x,  t)  . 


3x 


] =T(X, t) 


(2.20) 
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and  the  natural  boundary  condition  is 


[G(x)J(x)ii<5;-£>  ]  |x<=tfx(L,t) 


(2.21) 


For  free  torsional  vibration  r(x,t)  =  0,  thus 


G(x)J(x )  ]=Ip(*)  d 


(2.22) 


For  a  uniform  bar 


nr  d26(x,t)  _nrdzQ(x,t)  =0 


(2.23) 


Assume  a  solution  of  the  form 


0  (x,  t )  =4>(x)  cos  (Ci >t-a) 


(2.24) 


Substitute  equation  2.24  into  equation  2.23 


GJ  1-1^1 +<£  pi Jb(x)  =0 


(2.25) 


and  let 


P2  =  ffp<^ 

GJ 


(2.26) 


The  torsional  equation  of  motion  can  then  be  written  as 


-  *(x)  +P,»(x)  =0 


,2  r* 


(2.27) 


The  general  solution  is  of  the  form 


*  ( x )  s  in  ( Pxx )  +C2cos  ( pxx ) 


(2.28) 
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Numerical  Solutions 


The  position  dependent  equations  of  motion  for  x-y  plane  bending 
and  torsion  are  again 

d4V(f-~PSV(x)=0 

dx4  y 

■d2*<xUp^(x)=0 

dx 2 

with  general  solutions 

V (x) =A1sin( PyX) +A2cos ( PyX) +A3sinh ( PyX) +A4cosh(  PyX)  (2.31) 


(2.29) 


(2.30) 


4>(x)  =C1sin(  Pxx)  +C2cos  (  Pxx) 


The  boundary  conditions  at  the  clamped  end  are 


v(0,  t) =0 


8(0, t) =0 


av(0,t) -n 

i Jt 


At  x  =  L,  the  boundary  conditions  are 


wx=Idxx- 


dzQ(L,t) 

Tt*— 


dt‘ 


dv(x, t ) 


W 


(2.32) 

(2.33) 

(2.34) 

(2.35) 

(2.36) 

(2.37) 

(2.38) 


where  the  shear  force  and  moments  are 

After  applying  the  above  boundary  conditions  the  equation  of 
motion  for  x-y  plane  bending  reduces  to  2  equations  and  3  unknowns  while 
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(2.39) 


V  =EI  a3v(*, t) 
y  Z  8x3 

HX--GJ  a6(-X'-)  (2.40) 


tf2=£Iz 


d2v(x,t) 

dx 2 


(2.41) 


the  equation  of  motion  for  torsion  reduces  to  1  equation  and  2  unknowns. 
The  x-y  bending  plane  equation  of  motion  can  be  represented  in  matrix 
form  as 


^(Pyt)[Al  A4]T  =  [0] 


with 


*1  = 
*4= 


“As 

-a2 


F( PyL)  is  defined  as 


F(PyL)  = 


Fn 

F2'\  F22 


(2.42) 


(2.43) 


(2.44) 


where 

F11=-(pyL)2mp(sin(PyL)  -sinh(PyL) )  +£,pA(PyL)  (cos(  Pyl,) +cosh(  Pyl,)  )  (2.45) 
Fl2=(  PyI-)2nip(cos(  Pyi)  -cosh(  PyL) )  +LpA(  PyL)  (sin(  Pyl,)  -sinh( Py£) )  (2.46) 
F21  =  (  Py£ ) 3IZZ ( cos (PyL) -cosh (PyL)  )-L3pA(sin(PyL)  +sinh(PyL) )  (2.47) 

F22=(PyI,)3Izz(sin(PyL)+sinh(PyL)  )+L3pA(cos(PyL)+cosh(PyL) )  (2.48) 


For  equation  2.42  to  have  a  non-trivial  solution,  the  determinant  must 
be  zero.  The  natural  frequencies  for  bending  are  found  by  setting 
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detf ( PyL) =0 


(2.49) 


solving  for  PyL,  then  using 


(2.50) 


to  solve  for  the  natural  frequencies. 

After  applying  the  boundary  conditions  for  torsion  we  find  C2  =  0. 
The  equation  of  motion  can  be  written  as  a  transcendental  function 


tan( PXL) - 


P-TpL 

TO; 


=0 


(2.51) 


Of  PxL. 

The  natural  frequencies  are  found  from 


GJ 

PXP 


(2.52) 


Once  PyL  and  PXL  are  found,  the  constants  A.,,  A^,  and  C1  can  be  found 
from  equations  2.31  and  2.32. 

Equations  2.49  and  2.51  were  solved  numerically.  The  natural 
frequencies  of  the  beam  without  the  end  mass  were  also  determined. 

Table  2.3  contains  the  numerical  solutions  to  the  theoretical  equations 
of  motion. 
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Table  2.3  Numerical  Results 


Beam  without  End  Mass 

Beam  with  End  Mass 

XY  Bending 

M-l 

6.65 

2.17 

M-2 

41.72 

28.01 

M-3 

116.81 

76.36 

M-4 

228.91 

138.35 

M-5 

378.40 

240.29 

Torsion 

T-l 

391.89 

15.07 

T-2 

1176.00 

784.25 
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III.  Actuators,  Sensors,  Measurement  Channel,  and  Configuration 


The  original  Advanced  Beam  Experiment  conducted  by  Cristler  used 
two  pairs  of  linear  proof  mass  actuators  mounted  on  a  circular  disk 
attached  to  the  beam.  The  actuators,  mounted  in  orthogonal  axis, 
theoretically  allowed  for  simultaneous  control  of  the  beam  in  both 
bending  planes  and  torsion.  The  weight  of  the  actuators,  however, 
created  second  and  third  bending  modes  close  to  the  free  end  of  the 
beam.  In  addition,  the  actuators  could  not  achieve  full  force  output 
below  ^  Hz.  Since  the  first  bending  mode  in  both  the  x-y  plane  and  the 
x-z  plane  was  below  5  Hz,  effective  control  was  achieved  with  only  the 
first  torsion  mode. 

The  modified  ABE  conducted  by  Jacques  used  only  1  pair  of  linear 
proof  mass  actuators  to  control  first  mode  bending  in  the  x-y  plane  and 
the  first  torsion  mode.  To  control  second  and  third  mode  bending,  a 
structural  dynamics  shaker  was  added.  In  addition,  improved  accel¬ 
erometers  and  integration  circuits  were  used.  This  chapter  covers  the 
changes  implemented  in  the  modified  ABE  which  was  used  as  the  model  for 
this  work. 

Actuators 

In  the  original  ABE,  AFWAL  chose  linear  proof  mass  actuators  based 
(7:2.1)  on  a  design  by  TRW.  The  actuators  are  linear  dC  motors  which 
use  momentum  exchange  between  the  base  plate  and  the  moving  mass  to 
provide  a  control  force.  The  actuators  consist  of  a  linear  motor  coil 
mounted  on  two  support  brackets  connected  to  the  base  disk.  The  motor 
coil  is  driven  by  a  power  amplifier  circuit  which  transforms  a  voltage 
command  into  a  drive  current.  The  drive  current  is  limited  to  2  amps  by 
a  current  limiter  to  prevent  burning  out  the  motor  coils.  The  motor 
coils  drive  the  motor  magnets  which  are  contained  in  the  0.9258  kg 
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cylindrical  proof  mass.  The  proof  mass  travels  on  linear  bearings  which 
allow  for  +1.067E-02  meters  of  travel.  Below  5  Hz  the  maximum  force 
output  is  limited  by  the  distance  the  proof  mass  can  travel.  For  low 
frequencies,  the  peak  force  can  be  predicted  by 

^out  =  md^d  (3,1) 

where  md  is  the  mass  of  the  proof  mass  and  attached  accelerometer,  w  is 
the  frequency  in  rad/sec  of  the  signal  driving  the  motor,  and  d  is  the 
peak  displacement.  Figure  3.1  is  a  graph  of  the  theoretical  force 
output  for  each  actuator.  Above  5  Hz  the  maximum  force  output  of  an 
actuator  is  limited  to  about  8.9  N  by  the  motor  coil  capability  of  2 
amps.  Jacques  (7:18)  measured  the  actual  control  force  available  from 
the  actuators  to  control  mode  1  bending.  Table  3.1  gives  Jacques' 
measured  results  versus  predicted  results.  The  total  force  available  at 


Figure  3.1  Actuator  Theoretical  Force  Ouput 
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2  Hz  is  found  to  be  ±  2.625  N  which  is  the  sum  of  maximum  output  for 
both  actuators. 


Table  3.1  Actuator  Maximum  Force  Outputs  at  2  Hz  (7:18) 


Actuator 

Efficiency  (%) 

A 

1.606 

1.326 

82.5 

B 

1.606 

1.299 

80.9 

The  actuators  are  instrumented  with  Linear  Variable  Differential 
Transformers  (LVDT)  which  provide  feedback  of  the  relative  position 
between  the  proof  mass  and  the  motor  base.  The  proof  mass  is 
instrumented  with  an  Endevco  piezoresistive  accelerometer.  Figure  3.2 
is  a  schematic  of  the  linear  proof  mass  actuators  from  Christler 
(5:C.l  -  C. 4) .  Table  3.2  gives  a  description  of  each  component,  its 
mass,  and  its  dimensions.  Tables  3.3  and  3.4  list  the  center  of  gravity 


Table  3.2  Actuator  Mass  Model  Components  (5:C.l) 


Part 

Figure 

Mass 

Dimensions  (10"2m) 

Description 

(kg) 

Length 

Width 

Height 

Diameter 

1 

Rectangular 

Plate 

2.49E-2 

2.79 

0.76 

4.32 

- 

2 

Circular 

Cylinder 

3.49E-2 

11.43 

- 

- 

0.64 

3 

Rectangular 

Plate 

8.48E-2 

2.79 

0.76 

4.32 

- 

4 

Rectangular 

Plate 

3.49E-2 

2.79 

0.76 

5.40 

- 

5 

Circular 

Cylinder 

9.47E-2 

13.97 

- 

- 

1.91 

6 

Hollow 

Circular 

Cylinder 

51.36E-2 

5.72 

Inside  Diameter  = 
Outside  Diameter  = 

0.64 
=  5.40 
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and  mass  moments  of  inertia  for  each  component  of  actuators  A  and  B. 

The  specifications  for  the  linear  motors,  the  LVDT's,  and  the  Endevco 
accelerometers  are  listed  in  Appendix  B.  The  power  amplifier  circuit  is 
also  given. 


Table  3.3  Actuator  A  Mass  Model  (5:C.3) 


Part 

Center  of  Gravity  (1C)"2  m) 

Mass  Moment  of  Inertia 
(10'5  kg-m2) 

X 

y 

z 

xy 

1 

177.55 

-5.72 

-10.16 

0.17 

0.55 

0.40 

2 

176.78 

0.00 

|  1 

3.81 

0.018 

3.81 

3 

177.55 

5.72 

-10.16 

0.59 

1.87 

1.36 

4 

172.70 

5.72 

-10.16 

0.24 

1.07 

0.86 

5 

171.32 

6.20 

-10.16 

15.62 

0.43 

15.62 

6 

176.78 

0.00 

-10.16 

40.99 

33.13 

40.99 

Table  3.4  Actuator  B  Mass  Model  (5:C.4) 


Part 

Center  of  Gravity  (10‘2  m) 

Mass  Morrjpnt  of  Inertia 
(10’s  kg-m2) 

X 

y 

Z 

Jx 

Iz 

1 

177.55 

5.72 

10.16 

0.17 

0.55 

0.40 

2 

176.78 

0.00 

10.16 

3.81 

0.018 

3.81 

3 

177.55 

-5.72 

10.16 

0.59 

1.87 

1.36 

4 

172.70 

-5.72 

10.16 

0.24 

1.07 

0.86 

5 

171.32 

-6.20 

10.16 

15.62 

0.43 

15.62 

6 

176.78 

0.00 

10.16 

40.99 

33.13 

40.99 
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Initial  open  loop  testing  of  the  actuators  by  Cristler  found  them 
to  have  several  undesirable  characteristics .  These  characteristics 
(5:8-9)  are 

1.  The  low  frequency  behavior  was  non-linear  due  to  several 
factors.  The  bearing  friction  and  associated  hysteresis  caused 
drift  of  the  center  position.  Also,  the  limited  stroke  length 
would  not  allow  for  maximum  force  output  below  about  5  Hz.  For 
very  low  frequencies  (  <  1  Hz),  the  bearing  friction  often 
overcame  the  commanded  force  output  of  the  actuator. 

2.  The  zero  adjustment  of  the  power  amplifier  circuit 
required  continual  adjustment,  as  the  zero  position  would  shift  as 
a  function  of  both  frequency  and  amplitude. 

3.  The  open  loop  frequency  response  has  a  roll  off  and 
associated  phase  shift  in  the  vicinity  of  one  of  the  fundamental 
bending  mode  frequencies. 

4.  Non-linear  bracket  dynamics  appeared  in  the  region  of 
120-150  Hz. 

A  closed  loop  feedback  system  was  designed  to  overcome  these  open  loop 
characteristics.  The  details  of  this  design  are  completely  described  in 
Reference  7.  The  overall  goals  (7:2-4)  of  the  design  were: 

1.  Actuator  frequency  response  should  be  "flat”  over  the 
structure  control  bandwidth.  While  the  original  ABE  limited  the 
control  bandwidth  to  0  to  50  Hz,  the  modified  experiment  uses  an 
expanded  control  bandwidth  of  0  to  100  Hz. 

2.  The  actuator  proof  mass  should  maintain  an  inertia?,  position 
when  the  actuator  is  not  being  commanded. 

3.  The  proof  mass  should  maintain  its  centering  when  commanded  at 
different  frequencies  and  amplitudes. 

4.  The  actuator  proof  mass  should  remain  within  the  ±0.5  inch 
actuator  stroke  limits. 
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The  actuator  compensator  circuit  used  in  the  modified  ABE  was 
identical  to  the  circuit  used  in  the  original  experiment  except  for 
modifications  in  resistance  values.  These  modifications  (7:10)  were 
necessary  for  the  following  reasons: 

1.  The  characterization  of  the  open  loop  actuator  changed 
slightly.  The  low  frequency  pole  shifted  higher,  the  bearing 
friction  deteriorated  further,  and  the  actuators  were  sensitive  to 
environmental  changes. 

2.  The  modified  ABE  required  a  bandwidth  of  0-100  Hz  because  the 
third  mode  bending  frequency  of  60-70  Hz  was  included. 

3.  Different  feedback  accelerometers  were  used,  and  the 
difference  in  the  sensitivities  had  to  be  compensated  for. 

The  details  of  these  modifications  are  covered  in  detail  in  Reference  7. 
The  final  control  configuration  is  shown  in  Figure  3.3. 


Figure  3.3  Actuator  Compensation  Block  Diagram  (7:12) 
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Structural  Dynamics  Shaker 

The  third  actuator  was  added  to  control  second  and  third  mode 
bending.  The  actuator  used  was  an  Acoustic  Power  Systems  (APS)  Model 
113-LA  Structural  Dynamics  Shaker  driven  by  an  APS  Model  114  dual  mode 
power  amplifier.  This  shaker  was  chosen  because  of  its  long  stroke 
limit  and  excellent  response  in  the  control  bandwidth.  The 
exceptionally  flat  response  allowed  "open  loop"  operation,  thereby 
eliminating  the  feedback  problems  encountered  with  the  proof  mass 
actuators.  Since  the  structural  dynamics  shaker  was  physically  mounted 
in  the  laboratory,  the  original  intent  of  the  ABE  to  use  inertial 
controllers  mounted  to  the  beam  was  violated.  In  future  experiments, 
Jacques  suggested  (7:22)  that  the  structural  dynamics  shaker  should  be 
freely  suspended.  Since  the  goal  here  is  to  model  the  results  of 
Jacques'  work,  the  structural  dynamics  shaker  will  be  modeled  as  a  non- 
inertial  controller.  Details  of  the  open  loop  testing  are  in  Reference 
7.  Specification  sheets  for  the  shaker  and  the  amplifier  are  in 
Appendix  B. 

Sensors /Measurement  Channels 

One  of  the  goals  of  the  modified  ABE  was  to  improve  the  sensors 
and  integrator  circuit.  The  original  ABE  used  piezoelectric  accel¬ 
erometers  which  provided  relatively  poor  low  frequency  responses  and  the 
original  integration  circuit  was  susceptible  to  drift  caused  by 
integrating  dC  offsets  or  very  low  frequency  signals.  The  modified  ABE 
used  Endevco  Model  2262  piezoresistive  accelerometers  which  have  a  good 
low  frequency  response.  The  integrator  circuit  used  in  the  modified  ABE 
was  designed  by  WRDC/FIBG  to  have  negligible  phase  shift  and  a 
straight  -20  dB/decade  magnitude  slope  above  1  Hz.  Specifications  for 
the  accelerometers  and  the  integration  circuit  diagram  are  contained  in 
Appendix  B. 
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A  side  effect  of  signal  integration  was  attenuation  of  the  signal. 
Amplification  was  necessary  to  boost  the  analog-to-digital  input  signal 
up  to  the  minimum  range.  Two  amplification  stages  were  necessary.  The 
first  stage  was  placed  before  the  integrator  while  the  second  stage  was 
placed  after  the  integrator.  Gain  had  to  be  kept  low  in  the  first  stage 
to  avoid  clipping  the  Bignals  in  the  integrator.  The  remainder  of  the 
required  gain  was  provided  in  the  second  stage.  The  low  pass  filter  on 
both  amplifiers  was  set  to  1  KHz.  While  this  thesis  will  use  velocity 
feedback  as  shown  in  Figure  3.4  to  examine  the  integration  or 
structure/control  optimization,  the  control  law  will  be  implemented 
using  perfect  sensors. 
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Final  Configuration 

Figure  3.5  shows  the  modified  ABE  final  hardware  configuration. 

The  hardware  consisted  of  an  inverted  cantilevered  beam  with  base  plate, 
two  linear  proof  mass  actuators,  two  linear  proof  mass  actuator 
compensation  circuits,  a  structural  dynamics  shaker,  four 
accelerometers,  four  signal  conditioning/integrator  circuits,  a  PC-1000 
Array  processor  and  its  host  computer.  Figure  3.6  shows  the  system 
configuration  for  the  beam  with  base  plate,  the  actuators,  and  the 
structural  dynamics  shaker.  There  are  11  positions  located  7  inches 
apart.  The  shaker  is  located  at  position  4  (28  inches  from  the  fixed 
end)  and  provides  good  control  over  second  and  third  mode  bending.  The 
proof  mass  actuators,  shown  in  Figure  3.7,  are  mounted  parallel  to  the 
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Figure  3.6  ABE  Final  Configuration  (7:52) 


IV.  Finite  Element  Modeling  of  the  ABE 


The  finite  element  method  can  be  thought  of  as  a  version  of  the 
assumed-mode  method.  Instead  of  defining  a  function  which  describes  the 
deflected  shape  of  an  entire  structure,  the  finite  element  method 
defines  a  shape  function  f  over  a  small  region  in  such  a  way  so  as  to 
maintain  the  inter-element  continuity  of  the  i|»'s  for  the  entire 
structure.  A  brief  review  of  the  theory  from  Ref.  9  for  beams  in 
bending  and  torsion  is  presented  in  this  Chapter.  Ref.  10  contains  a 
more  detailed  discussion  of  finite  element  theory.  Also  included  in 
this  Chapter  is  the  development  of  the  ABE  finite  element  model  using 
MSC/NASTRAN. 

Theory 

If  we  assume  Bernoulli-Euler  beam  theory,  the  transverse 
displacements  for  a  beam  in  bending  can  be  approximated  by  (9:385-387) 

4 

v(x,t)=£  *j(x)v(t)  (4.1) 

i«1 

where  the  shape  functions  *  satisfy  the  boundary  conditions 


*1(0)  =  1, 

Vi  (0)  =  *,(£)  =  V^L)  =  0 

(4.2) 

v'z(O)  =  1 

*2(0)  =  *2  (L)  =  */2(i)  =  0 

(4.3) 

*5  (L)  =  1/ 

*j(0)  =  VjfO)  =  1|/3(L)  =  0 

(4.4) 

V«<£)  = 

*4(0)  =  */4(0)  =  *4(L)  =  0 

(4.5) 

For  a  beam  in  bending  the  general  solution  for  static  deflection  v(x)  is 
the  cubic  polynomial 
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’<*>  ■  *  Mi)  *  Mi)!  *  Hi)3 


(4.6) 


After  substituting  the  given  boundary  conditions  into  equation  4.6,  the 
beam  element  shape  functions  are  found  to  be 


*’  ■ 1  -  iif  *  2(if 


(4.7) 


*  -  *  -  MU  *  izf 


(4.8) 


*  ■  Hf  -  W 


(4.9) 


*  ■  -iif  *  ^)3 


(4.10) 


By  applying  equation  4.6  to  the  kinetic  and  potential  energy 
expressions,  terms  for  the  stiffness  and  mass  are 


k. 


(4.11) 


=JpA\|rii|rjC 


(4.12) 


Substituting  the  shape  functions  into  the  expressions  for  the  stiffness 
and  mass  results  in  the  element  stiffness  and  mass  matrix 


U3J  - 


12 

6  L 

-12 

6  L 

6  L 

4i2 

-6L 

2  L2 

-12 

-6  L 

12 

-6  L 

6  L 

2  L2 

-6  L 

4  L2 

(4.13) 


4-2 


(4.14) 


[ml 


156 

22L 

54 

-13  L 

22L 

4  L2 

13L 

-3  L2 

54 

1 3  L 

136 

-22L 

-13L 

-3  L2 

-22  L 

4  L2 

The  torsional  stiffness  and  mass  matrix  for  a  beam  element  can  be 
found  in  a  similar  manner  (9:388-389).  The  rotation  along  the  element 
is  given  by 

0(x,t)  =4r1(x)01(t)  +  T^(x)02(t)  (4.15) 

subject  to  boundary  conditions 

»M0)=1,  *i(£)=0  (4.16) 


ife(0)=0,  *2(L)= 1  (4.17) 

The  shape  functions  are  given  by 

t,(x)-l-*  (4.18) 

Li 

and 

i|r2(x)=J  (4.19) 

Substituting  equation  4.15  into  the  torsional  expressions  for  potential 
and  kinetic  energy,  the  expressions  for  stiffness  and  mass  become 


l 


kjj  =  [GJtj^jdx 


(4.20) 


r 

=  P 


mij  ~  plpWjd* 


(4.21) 
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Inserting  the  shape  functions  from  equations  4.18  and  4.19  into 
equations  4.20  and  4.21  gives  the  stiffness  and  mass  matrix  for  a  beam 
element  in  torsion. 


(4.22) 


[m] 


(4.23) 


For  an  undamped  system,  the  element  matrices  can  be  formed  into 
partitioned  system  mass  and  stiffness  matrices  (9:406-407).  The 
matrices  are  partitioned  into  active  and  constrained  degrees-of-freedom 
giving  the  equation  of  motion 


■ 

/  .. 

/ 

f 

^aa 

^ac 

+ 

^aa 

■^ac 

= 

Pa 

^ca 

Mcc 

Pc. 

■^ca 

•Kcc 

Uc 

pc. 

for  the  system.  For  Uc  =  0,  equation  4.24  can  be  reduced  to 


MaaUa+KasUa=Pa 


(4.25) 


plus  an  additional  equation  of  reaction  constraints.  Since  Mga  and  Kga 
are  sufficient  to  solve  for  the  active  displacement  vector  Ug,  the 
second  equation  is  not  required. 

Assume  a  solution  of  the  form 

(U)  =(^r)coswt  (4.26) 

where  (i|i)  is  a  vector  of  real  numbers  and  cos  (cot)  is  a  scalar 
multiplier.  Substituting  equation  4.26  into  equation  4.25  gives 

[K-<Jm]  (\}f) coscot  =0  (4.27) 

For  this  equation  to  be  valid  for  all  time. 
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[K-(£m]  (t|t)  =0 


(4.28) 


The  solution  to  this  equation  gives  the  eigenvalues  e,  and  their 
associated  eigenvectors. 

MSC/NASTRAN  66A  Solution  Solver 

MSC/NASTRAN  solution  sequence  "Real  Eigenvalue  Analysis"  (SOL  3) 
was  used  for  the  preliminary  modeling  of  the  ABE.  "Modal  Transient 
Response"  (SOL  112)  was  used  in  the  optimization  sequence  so  that 
initial  conditions  could  be  obtained.  In  all  cases,  the  eigenvalue 
extraction  method  used  was  the  Modified  Givens  Method  with  mass 
normalization.  Unlike  the  Givens  method,  the  Modified  Givens  Method 
does  not  require  that  the  mass  matrix  be  nonsingular.  Instead  of 
performing  Cholesky  decomposition  on  the  mass  matrix,  a  positive 
definite  matrix  [K  +  A.MJ  is  formed  where  X  is  selected  by  the  program  to 
optimize  the  reliability  and  accuracy  of  the  eigenvalue  extraction.  The 
Givens  transformation  method  is  used  to  convert  the  program  defined  [J] 
matrix  to  tridiagonal  form.  A  modified  Q-R  algorithm  is  used  to  extract 
the  eigenvalues  from  the  tridiagonal  matrix.  The  complete  eigenvalue 
extraction  procedure  is  contained  in  Reference  (11:4.2-4  -  4.2-8). 

Initial  Beam  Model 

SDRC-IDEAS  was  used  (12)  as  the  pre/post-processor  for  the 
MSC/NASTRAN  data  file.  The  ABE  finite  element  model  consists  of  11  grid 
points,  or  nodes,  which  are  connected  by  NASTRAN  CBAR  elements.  The 
CBAR  element  specifies  the  beam  connecting  nodes  and  references  a 
physical  and  a  material  property  table  for  the  element.  The  SPC  entry 
is  used  to  specify  the  physical  degrees  of  freedom  for  each  grid  point. 
For  the  ABE  model  in  this  thesis,  the  fixed  end  has  no  degrees  of 
freedom.  All  other  nodes  are  free  to  translate  in  the  y-direction  and 
rotate  about  the  "x"  and  "y"  axis. 
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The  first  run  model  consisted  of  a  10  element  beam  without  the  end 


disk  or  actuators.  Results  were  computed  in  both  SDRC-IDEAS  and 
MSC/NASTRAN  and  compared  to  the  numerical  calculations  performed  in 
Chapter  II.  Table  4.1  gives  the  natural  frequencies  calculated  by  using 
both  programs.  For  the  simple  beam  without  the  end  disk,  MSC/NASTRAN 
quickly  losses  accuracy.  This  loss  of  accuracy  for  a  simple 
cantilevered  beam  using  MSC/NASTRAN  is  also  shown  in  Reference  (13:7.9). 
No  reason  is  given  in  the  reference  for  such  a  large  error  in  the  third 
and  higher  bending  modes.  The  error  is  reduced  considerably  once 
the  base  disk  is  attached. 


Table  n . 1  Beam  Without  End  Disk 


Mode 

Numerical 

MSC/NASTRAN 

Error 

SDRC-IDEAS 

.1 

Error 

Hz 

Hz 

% 

Hz 

% 

XY-1 

6.65 

6.61 

-  0.60 

6.65 

0.00 

XY-2 

41.72 

40.88 

-  2.01 

41.65 

-  0.17 

XY-3 

116.81 

112.81 

-  3.42 

116.44 

-  0.32 

T— 1 

391.89 

353.59 

-  9.77 

392.40 

+  0.13 

The  second  step  in  the  model  construction  was  to  construct  a  rigid 
body  disk  in  SDRC-IDEAS  and  attach  it  to  the  beam.  A  concentrated  mass 
model  was  also  constructed  using  the  inertia  properties  of  the  disk  and 
as  expected  the  dynamic  results  were  identical.  As  a  result,  the 
concentrated  mass  model  was  used.  Table  4.2  gives  the  natural 
frequencies  of  the  beam  with  the  end  disk  attached  using  both  SDRC-IDEAS 
and  MSC/NASTRAN  and  compares  their  results  to  the  numerical  results 
obtained  in  Chapter  II. 
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Table  4.2  Beam  With  End  Disk 


Mode 

Numerical 

MSC/NASTRAN 

Error 

SDRC/ IDEAS 

Error 

Hz 

Hz 

Hz 

% 

XY-1 

2.17 

2.17 

0.00 

2.17 

0.00 

XY-2 

28.01 

27.92 

0.32 

27.99 

-  0.07 

XY-3 

76.36 

76.17 

-  0.25 

76.25 

-  0.14 

T-l 

15.07 

15.05 

-  0.13 

15.05 

-  0.13 

As  mentioned  earlier,  the  disk  provides  the  primary  component  of 
rotatory  inertia.  From  Table  4.2,  it  is  apparent  that  the  torsion  mode 
is  now  well  within  the  0-100  Hz  bandwidth  of  the  modified  ABE.  Figures 
4.1  -  4.4  show  the  mode  shape  comparisons  for  the  first  three  bending 
modes  and  torsion.  The  top  figure  on  each  page  is  the  mode  shape  of  the 
beam  without  the  end  disk.  The  bottom  figure  is  the  mode  shape  for  the 
beam  with  the  end  disk.  Notice  the  mode  shape  for  mode  2  and  3.  The 
weight  of  the  disk  is  causing  the  beam  to  act  almost  like  a  -inned  beam 
at  the  free  end.  This  was  noticed  in  earlier  work  (7:4)  on  the  ABE  and 
was  one  of  the  reasons  for  poor  control  results  for  these  modes.  The 
modified  ABE  uses  the  structural  dynamics  shaker  to  overcome  this 
problem. 

The  third  step  was  to  add  the  actuator  non-moving  mass  components 
to  the  beam.  Tables  3.2,  3.3,  and  3.4  contain  the  mass  and  inertia 
properties  of  the  actuator  components.  Once  again,  these  were  added  to 
the  finite  element  model  as  concentrated  masses  with  inertial 
properties.  The  results  of  this  step  are  listed  in  Table  4.4.  This 
step  is  referred  to  by  Jacques  as  the  "clean  configuration." 

The  final  configuration  includes  the  actuator  mass.  This 
configuration  was  not  developed  in  SDRC-IDEAS.  The  natural  frequencies 
for  this  configuration  are  presented  in  a  later  section  in  Table  4.4. 
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Figure  4.1  Mode  1  Bending  For  Beam  and  Beam  With  Disk 
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Figure  4.3  Mode  3  Bending  For  Beam  And  Beam  With  Disk 
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Figure  4.4  Torsion  Mode  For  Beam  And  Beam  With  Disk 
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Table  4.3  Clean  Configuration 


Model  Verification 

Even  though  the  ABE  is  a  very  simple  structure,  the  finite  element 
model  created  in  MSC/NASTRAN  still  contains  modeling  errors.  The 
equations  of  motion  from  the  finite  element  model  only  represent  a 
mathematical  model  of  the  real  structure.  The  goal  in  modeling  is  to 
match  as  closely  as  possible  the  behavior  of  the  physical  structure  in 
the  realm  of  interest.  The  realm  of  interest  for  the  ABE  is  the  dynamic 
behavior  of  the  first  three  bending  modes  and  the  first  torsion  mode. 

Convergence.  To  determine  if  10  elements  were  sufficient,  10  beam 
models  were  created  with  different  numbers  of  elements  to  test  for 
convergence  of  the  natural  frequencies.  Table  4.4  contains  the  results 
of  this  test. 


Table  4.4  Convergence  Test  Results 


Beam  - 
Number  of 
Elements 

XY-1 

XY-2 

XY-3 

T-l 

B-2 

2.174 

28.182 

79.088 

15.045 

B-3 

2.174 

28.029 

76.848 

15.045 

B-4 

2.174 

28.000 

76.447 

15.045 

B-5 

2.174 

27.992 

76.333 

15.045 

B-6 

2.174 

27.989 

76.290 

15.045 

B-7 

2.174 

27.988 

76.271 

15.045 

B-8 

2.174 

27.987 

76.262 

15.045 
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Beam  - 
Number  of 
Elements 

XY-1 

XY-2 

XY-3 

T-l 

B-9 

2.174 

27.987 

76.256 

15.045 

B-10 

2.174 

27.986 

76.253 

15.045 

B-35 

2.174 

27.986 

76.246 

15.045 

Numerical 

2.174 

28.010 

76.360 

15.072 

Error  B-10 
<%) 

0.00 

-0.09 

-0.14 

-0.18 

The  finite  element  model  converges  from  above  to  a  solution  below 
the  numerical  results.  For  the  10-element  beam,  the  accuracy  is  quite 
good  for  the  first  three  bending  modes  and  the  first  torsion  mode. 

Effects  of  Shear  Deformation.  The  theoretical  model  was  based  on 
Bernoulli-Euler  beam  assumptions.  These  assumptions  neglected  shear 
deformation  and  rotatory  inertia.  From  (9:204),  the  theoretical 
equation  of  motion  for  a  short,  stubby  beam  (Timoshenko  Beam  Theory) 


effects 

is 

/ 

ei^L  - 

- 

pi  sf 

dx4 

„  dt 2 

.  yi  r 

dx2dt 2 

_  ( 

El  dt 


K  GA 


dx‘ 


p  -  pA 


d2v 
dt 2 


pj  dr 


KGA 


dt‘ 


p  -  pA 


ddv 
dt 2) 


(4.29) 


=  0 


The  first  term  in  brackets  is  from  Bernoulli-Euler  beam  theory.  The 
principal  rotatory  inertia  is  in  the  second  term  while  the  principal 
shear  deformation  is  in  the  third  term.  The  fourth  term  includes 
combined  rotatory  inertia  and  shear  deformation.  SDRC-IDEAS  accounts 
for  the  effects  of  shear  by  the  definition  of  a  shear  area  ratio  term. 
Ihe  shear  area  ratios  are  defined  in  SDRC-IDEAS  to  be  (12:2.18) 
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(4.30) 


a/2  02 

Sr„=^_  f  Sidy 


' ry 


ly 


l 


2  1  a 


b/2  q2 


where,  the  Q's,  shown  in  figure  4.5,  are 

a/2  b/2 

Qy=  f  f  zdzdy 
-1/2  0 

b/2  a/2 

02=  f  f  ydydz 
-b/2  0 


Figure  4.5  Q  Area  For  Shear  Area  Ratio  Calculation 


(4.31) 


(4.32) 


(4.33) 
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Results  from  the  10-element  beam  model  were  computed  with  and  without 
the  effects  of  shear.  Table  4.5  presents  the  results  of  these 
computations.  The  data  suggests  that  the  Bernoulli-Euler  assumption  was 
a  good  assumption.  For  the  final  model,  however,  the  effects  of  shear 
and  rotatory  inertia  will  be  included. 


Table  4.5  Effects  of  Shear 


Mode 

B-1Q  W/O  Shear  (Hz) 

B-10  With  Shear  (Hz) 

Comparison  (%) 

XY-1 

2.174 

2.174 

0.00 

XY-2 

27.993 

27.986 

-0.03 

XY-3 

76.280 

7  6 . 246 

-0.04 

T-l 

15.045 

15.045 

0.00 

Comparison  With  Modal  Test  Data.  The  final  "test"  of  the 
mathematical  finite  element  model  is  to  compare  it  to  experimental  data. 
In  structural  dynamics  modeling,  finite  element  analysis  and  modal 
analysis  take  different  approaches  to  obtain  a  system  model  (14;86). 
Modal  analysis  models  are  derived  experimentally  while  finite  element 
models  are  derived  mathematically.  This  difference  makes  modal  analysis 
a  good  check  of  the  finite  element  model.  Jacques  performed  modal 
analysis  on  the  ABE  to  determine  its  response  to  vibration.  Table  4.6 
gives  the  comparison  between  Jacques  modal  analysis  model  of  the  final 
configuration  (7:39)  and  the  MSC/NASTRAN  finite  element  model  with 
actuators  and  disk.  One  reason  for  the  large  difference  in  second  and 
third  mode  bending  may  be  because  displacements  in  the  x-direction  were 
constrained  in  MSC/NASTRAN.  Another  reason  for  the  difference  is 
because  of  the  exclusion  of  the  actuator  dynamics  from  the  finite 
element  model.  The  section  on  damping  will  address  the  actuator 
dynamics . 
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Modal  Transformation 


The  finite  element  model  equation  of  motion  which  includes  the 
mass,  stiffness,  and  damping  matrices  as  well  as  the  control  inputs  is 

Table  4.6  Final  Configuration  Model  Comparison 


Mode 

Jacques  Modal  Test  (Hz) 

MSC/NASTRAN  (Hz) 

Error  (%) 

XY-1 

1.97 

1.93 

-  2.03 

XY-2 

23.43 

27.74 

+  18.40 

XY-3 

62.92 

73.62 

+  14.53 

T— 1 

13.35 

13.07 

-  2.10 

[m](q)  +  [c](4)  +  [k](q)  =  [£>]  (u)  (4.34) 

where  q  is  a  vector  of  generalized  coordinates,  [b]  is  the  control  input 
distribution  matrix,  and  u  is  a  vector  of  control  forces.  Assuming  that 
only  the  first  two  bending  modes  and  the  first  torsion  mode  are  present 
in  the  response,  then  q  is  a  3x1  vector  and  the  generalized  coordinates 
are  defined  as 

ql  =  displacement  of  node  11. 
q2  =  rotation  of  node  11. 
q3  =  displacement  of  node  5. 

In  terms  of  the  physical  coordinates,  the  generalized  coordinates  are 

=  -j  (  ~/l1a  +  y  11b>  (4.35) 


?  2  = 


+  y  11b> 


(4.36) 


<?  3  =  ys 


(4.37) 


In  terms  of  the  generalized  coordinates,  the  measured  coordinates  are 
then 
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yila  =  -?1  +  $2 


(4.38) 


yilb  =  4l  +  72  (4.39) 

y5  =  g3  (4.40) 

If  we  assume  the  modes  are  uncoupled  and  assume  modal  damping,  equation 
4.34  can  be  decoupled  by  using  modal  coordinates  t)  and  the  modal  matrix 
of  right  eigenvectors  *  to  define  a  transformation  such  that 

[g]  =  *  [til  (4.41) 

where  the  eigenvectors  have  been  normalized  such  that 

*T[m]®  =  [tf]  =  [I]  (4.42) 

Equation  4.34  can  then  be  written  as 

(ii)  +  (201(11)  +  [  cc^]  (n)  =  ®T[i>]  (u)  (4.43) 

The  matrix  [2£a>]  is  a  diagonal  matrix  where  the  ith  element  represents 
damping  for  the  ith  mode.  The  matrix  [o2]  is  a  diagonal  matrix  where 
the  diagonal  entries  are  the  eigenvalues  for  equation  4.43. 

Damping 

Damping  is  inherent  in  any  structure,  but  predicting  how  much 
damping  is  in  a  structure  is  difficult.  Modal  testing  provides  a  good 
approach  to  help  understand  the  effects  damping  has  on  a  structure. 

From  modal  analysis  data,  Jacques  estimated  how  much  damping  was  present 
in  the  ABE.  Table  4.6  gives  the  experimentally  determined  damping 
values  for  both  the  "clean  configuration"  and  the  final  configuration 
with  actuators  attached.  The  viscous  damping  factor  C  is  the  ratio 
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c / ccr  where  ccr  is  the  critical  damping  factor.  For  large  space 
structures,  the  damping  factor  is  typically  very  small. 

By  including  the  actuator  damping  with  the  beam  damping,  the 
stability  effects  of  the  actuators  cannot  be  determined.  Inman  (15:508- 


Table  4.6  Modal  Analysis  Damping  Values  (7:39,41) 


Mode 


XY-1 


XY-2 


XY-3 


T-l 


Clean  Configuration  (Q 


0.00129 


0.00375 


0.00086 


Not  Determined 


Final  Configuration  ({) 


0.0640 


0.0236 


0.0121 


0.00839 


Actuator  Mass 


Norunoving  Mass  of 
the  Actuator 

Structure 


ka  '  Actuator  Stiffness 
ca  -  Actuator  Damping 
fa  -  Actuator  Force 


Figure  4.6  Simple  Model  of  Actuator  Dynamics  (15:508) 
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513)  shows  the  effects  of  actuator  dynamics  on  the  overall  system 
stability.  Figure  4.6  shows  a  simple  model  of  the  actuator  dynamics. 

If  actuator  dynamics  are  neglected,  the  mathematical  "answer"  to  how 
much  velocity  feedback  gain  should  be  added  to  achieve  a  desired 
response  is  as  much  as  is  necessary.  Inman  presents  examples  which  show 
how  the  system  can  go  unstable  as  the  gain  is  increased.  When  the 
actuator  dynamics  are  introduced  in  the  single-degree-of-freedom  model 
shown  in  Figure  4.6,  the  [C]  and  [K]  matrices  in  the  equations  of  motion 
become  coupled.  If  velocity  feedback  is  used  like  it  is  in  the  ABE,  the 
[C]  matrix  becomes  coupled  and  unsymmetric.  This  limits  the  amount  of 
gain  that  can  be  added  before  the  actuators  cause  the  system  to 
destabilize. 

Since  this  thesis  is  modeling  the  modified  ABE,  actuator  damping 
will  be  treated  as  viscous  modal  damping.  Proportional,  or  Rayleigh, 
damping  can  be  used  to  model  the  damping  where 

(C]  =  af  K]  +  P[W]  (4.44) 

Matrix  [C]  is  an  orthogonal  damping  matrix  (10:377)  because  it  permits 
modes  to  be  uncoupled  by  the  eigenvectors  associated  with  the  undamped 
eigenvalue  problem.  The  damping  ratio  is  related  to  the  stiffness  and 
mass  proportional  damping  constants  a  and  P  by 


a  and  P  can  be  found  in  terms  of  u  and  (  by 

a  =  2((2<o2-C1a>1)/(a|-<^)  (4.46) 

P  =  2g>]G$((']Ci£-(2(‘>i)/(c‘^-c4)  (4.47) 

Using  data  in  Table  4.6,  a  =  2.2072E-04  and  P  =  1.5219.  The  graph  for 
0  <  <■>  <  400  rad/sec  is  shown  in  Figure  4.7  where  the  Rayleigh  damping 
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curve  is  the  sum  of  the  mass  proportional  damping  curve  and  the 
stiffness  proportional  damping  curve. 


Figure  4.7  Rayleigh  (Proportional)  Damping  Model 


Initial  Conditions 

Structures  are  normally  designed  to  withstand  some  given  design 
conditions.  When  the  design  conditions  are  applied,  dynamical  systems 
will  exhibit  a  response.  This  response  is  needed  as  initial  conditions 
for  each  structural  iteration  in  the  structural/control  optimization. 
Miller,  Vipperla  and  Venkayya  (3:928)  demonstrated  that  it  is  not 
adequate  to  specify  the  initial  response  conditions  for  a  structure 
being  optimized  and  then  hold  these  constant  for  the  iteration  process. 
As  the  structural  property  is  varied  (cross-sectional  area)  in  the 
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optimization  process,  each  iterated  structure  will  exhibit  different 
initial  conditions.  Therefore,  it  is  necessary  to  specify  the  design 


conditions,  such  as  the  applied  load,  and  then  solve  for  the  initial 
conditions  for  each  iteration  of  the  structure. 

An  impact  load  was  included  in  the  MSC/NASTRAN  finite  element 
solution  process.  To  obtain  the  initial  conditions,  the  TLOAD  data 
entry  is  used  to  specify  a  load  table  TABLED 1  and  a  scale  factor  entry 
DAREA.  The  equation  for  the  force  input  in  MSC/NASTRAN  is 

P(t)  =  A[E( t -r ) ]  (4.48) 

where  the  time  delay  r=0.  The  scale  factor  A  was  selected  to  be  1000  N. 
The  load  table  is  shown  graphically  in  Figure  4.7.  The  table  is  shown 
graphically  to  illustrate  how  MSC/NASTRAN  reads  input  tables.  If  the 
three  points  at  t=0.001,  t=0.002,  and  t=0.003  where  the  only  points 


Figure  4.8  MSC/NASTRAN  Force  Area  Chart 


4-21 


specified,  an  infinite  load  would  be  applied.  This  is  because 
MSC/NASTRAN  extrapolates  tables.  Since  no  value  is  specified  after 
t=0.003,  the  slope  from  t=0.002  to  t=0.0Q3  would  be  continued  after 
t=0.003.  To  prevent  an  infinite  load,  values  of  zero  must  be  specified 
on  both  sides  of  the  load.  The  initial  conditions  for  the 
structure/control  optimization  were  chosen  at  t=0.003.  This  has  the 
effect  of  "turning  on"  the  control  system  at  t=0.003.  The  ABE  undamped 
and  damped  response  to  the  initial  conditions  are  shown  in  Figures  4.9 
and  4.10. 
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Figure  4.9  Response  of  Undamped  ABE  to  Impulse  Force 
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Figure  4.10  Response  of  Damped  ABE  to  Impulse  Force 


V.  Control  Theory 


The  state-space  formulation  is  detailed  in  this  Chapter.  LQR 
theory  and  the  selection  of  the  Q  and  R  weighing  matrices  are  presented. 


State-Space  Formulation 

After  normalizing  the  eigenvectors  with  respect  to  the  mass 
matrix,  the  linear  differential  equations  of  motion  from  equation  4.43 
can  be  written  in  state-space  form.  The  standard  state-space  format  for 
the  state  equation  is 

[x] =A[x] +B[u] +[w]  (5.1) 

where  the  state  vector  is 


[ x  ]  T  =  [  >  *12'  ^3  >  ^1 '  *12 '  ^3 1 

and  the  input  vector  is 

[  u]T  =  (f’A,FB,Fc] 

The  plant  matrix  A  and  the  input  distribution  matrix  B  are 

0  I 
A=  , 

-or  -2$o) 


(5.2) 


(5.3) 


(5.4) 


4>T£> 


(5.5) 


and  [w]  is  a  process  noise  or  input  disturbance  vector. 
The  output  equation  can  be  written  as 


[y)=C[x]+D[u]+[v] 


(5.6) 
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The  output  distribution  matrix  C  is 


C=[Cp*|cv»] 


(5.7) 


where  the  partitions  Cp  and  Cy  are  the  position  and  velocity  sensor 
distributions  respectively.  The  matrix  D  allows  the  output  variables  to 
be  controlled  by  [u]  while  [v]  is  a  measurement  noise  term.  The 
schematic  of  the  state-space  model  is  shown  in  Figure  5.1. 

This  thesis  will  assume  that  the  plant  disturbance  or  process 
noise  term  [w]  is  negligible.  In  addition,  the  sensors  are  assumed  to 
be  perfect  and  there  is  no  control  over  the  output  vector.  The  implies 
that  [v]=0  and  D=0.  The  input  distribution  matrix  [b]  is  defined  to  be 


-110 
110 
0  0  1 


Since  position  feedback  is  not  used,  Cp=0.  Cv  is  given  by 


(5.8) 


-110 
110 
0  0  1 


(5.9) 


thus,  for  collocated  sensors  and  controllers,  B=CT.  The  final  state- 


space  formulation  is 


[  x  ]  =A  [  x  ]  +B  [  u  ] 


(5.10) 


[y]=C[x] 


(5.11) 


Linear  Quadratic  Regulator  Theory 

Ref.  16  covers  LQR  theory  (16:6.1-6.9)  in  detail.  A  brief  review 
of  the  theory  is  presented. 

If  full  state  feedback  is  used,  the  control  input  can  be  defined 
as  a  linear  combination  of  the  state  variables 
that  make  the  closed-loop  system 
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Disturbances 


Commands  Noise 


Figure  5.1  Control  Schematic  (1:519) 


u  ( t )  =  -Kcx  ( t ) 


(5.12) 


x(t)=Ax(t)+fl[-Kcx(t)  ] 


(5.13) 


X ( t )  =  [A -BKC ] x ( t ) 


(5.14) 


asymptotically  stable.  The  control  gain  Kc  can  be  found  by  minimizing 
the  performance  index 


(xT(t)Ocx(t)+uT(t)Acu(t) ] 


dt 


(5.15) 


where  Qc  and  Rc  are  symmetric,  real,  positive  definite  matrices  chosen 
by  the  design  engineer  according  to  the  importance  of  the  state  and  of 
the  controls.  From  (16:6.1-6.2),  this  can  be  interpreted  as  "we  wish  to 
find  a  control  law  u(t)  such  that  the  integral-squared-error  of  the 
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deviation  of  the  state  trajectories  from  their  nominal  are  kept  small 
without  using  a  great  deal  of  control  energy."  From  optimal  control 
theory,  the  gain  matrix  Kc  in  equation  5.12  which  minimizes  equation 
5.14  is  given  by 

Kc=R~‘iB1S  (5.16) 

where  S  is  the  solution  to  the  algebraic  Riccati  equation 

ATS+SA-SBRc'1STS+Oc=0  (5.17) 

Since  position  feedback  is  not  used  for  this  work,  however,  the  gain 
matrix  Kc  is  not  guaranteed  to  be  optimal.  This  is  the  basis  of  a 
Linear  Quadratic  Regulator. 

Generally,  the  state  variables  cannot  be  measured  directly.  An 

estimator  is  required  to  reconstruct  the  states  from  the  sensor  outputs. 
The  estimator  has  the  form 

[  ■£  1  =A  [  x  ]  +B  [  u  ]  +K0  [y-y]  (5.18) 

[y]=C[x]  (5.19) 

where  [x]  and  [y]  are  the  estimated  state  vector  and  the  estimated 
output  vector.  The  observer  gain  KQ  must  be  chosen  so  that  the  error 
defined  by 

[e]=[x]-[x]  (5.20) 

is  stable.  Since  the  state  vector  is  not  available,  the  control  input 
is  also  baeed  on  the  estimated  state  vector 

[ u ] = -K0  [  x  )  (5.21) 

Kq  is  selected  to  minimize  the  estimator  cost  function 
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J0=|[xT(t)00x(t)+yT(t)^oy(t)]dt  (5.22) 

where  Q0  represents  the  process  noise  covariance  and  RQ  represents  the 
measurement  noise  covariance.  This  is  the  basis  of  a  Kalman  filter. 

When  a  Kalman  filter  is  added  to  a  Linear  Quadratic  Regulator,  the 
result  is  a  Linear  Quadratic  Gaussian  compensator. 


Force  Actuators  Flexible  System  Sensors 


N  -  Controlled  Modes 
R  -  Residual  Modes 


Figure  5.2  Control  and  Observation  Spillover  (17:425) 

If  residual  modes  are  present,  they  will  result  in  spillover. 
Figure  5.2  shows  a  schematic  of  spillover.  If  the  sensor  outputs  are 
contaminated  by  the  residual  modes  through  CR,  the  effect  is  known  as 
observation  spillover.  If  the  control  actuators  are  excited  by  the 
residual  modes  through  BR,  the  effect  is  known  as  control  spillover. 

The  effects  of  spillover  can  be  detrimental.  Calico  and  Janiszewski 
(18)  showed  that  eliminating  either  observation  or  control  spillover  was 
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sufficient  to  ensure  stability  of  the  suppressed  modes.  Jacques  (7) 
attempted  to  demonstrate  modal  suppression  on  the  ABE  by  eliminating 
observation  spillover.  The  results  achieved  while  trying  to  suppress 
bending  mode  2  were  partially  successful.  The  reason  for  the  partial 
success  was  attributed  to  the  estimator  model.  Balas  (17:421)  presents 
two  estimator  designs.  The  first  is  based  on  the  Kalman  filter.  The 
Kalman  filter  is  used  by  Balas  if  the  system  signal-to-noise  ratios  are 
not  very  high.  If  the  system  signal-to-noise  ratios  are  sufficiently 
high,  the  system  is  treated  as  deterministic  and  the  estimator  used  is  a 
Luenberger  observer.  The  reason  for  this  is  the  estimator  error  term  in 
equation  5.19  will  contain  additive  noise  terms  with  the  Kalman  filter 
but  these  terms  are  absent  with  the  Luenberger  observer.  Since  system 
noise  was  attributed  in  part  to  poor  estimator  performance,  additional 
work  on  the  ABE  should  examine  the  effect  of  estimator  choice  on  the 
control  model. 

For  the  integrated  structure/control  optimization  procedure 
implemented  in  this  thesis,  it  is  assumed  that  the  state  is  completely 
observable.  LQR  theory  will  be  used  as  the  optimal  control  law  without 
the  introduction  of  a  Kalman  filter.  MATLAB  uses  eigenvector 
decomposition  of  an  associated  Hamiltonian  matrix  to  solve  the  LQR 
problem  (19:CR-38,CR-39) .  For  convenience,  the  "c"  and  "o"  subscript 
will  dropped.  Parameters  used  for  the  remainder  of  this  work  refer  to 
the  controller. 

Choice  of  Weighing  Matrices  0  and  R 

In  the  past,  an  optimized  structure  was  given  to  a  control 
engineer  to  design  an  optimal  control  system.  The  engineer  selected  Q 
and  R  subject  to  the  importance  of  the  state  or  the  control  system  so 
that  a  desired  dynamic  response  was  achieved.  Venkayya  and  Tischler 
(20:433-434)  have  proposed  that  the  selection  of  Q  and  R  should  not  be 
arbitrary,  but  that  their  selection  should  ensure  that  the  cost  function 
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J  provides  a  measure  of  some  appropriate  physical  quantity,  such  as 
total  system  energy. 

From  Chapter  IV,  the  finite  element  equation  of  motion  is 

mq(t) +cg(t) +kg(t) =bu(t)  (5.23) 

If  J  is  defined  as 


J  =J  [  0mgTmg+0|cgTkg+0RuTjbTk~1jbu  ]dt  (5.24) 

0 

for  positive  scaling  parameters  0m,  0k,  and  0R,  then  J  is  the  absolute 
sum  of  the  kinetic,  strain,  and  potential  energies.  Total  system  energy 
is  thus  minimized  to  satisfy  the  LQR  control  law.  Considering  the  modal 
coordinate  transformation  4.41  and  the  mass  normalization  4.42,  Q  and  R 
can  be  written  as 


OnrT  0 

0  0k[«?) 


(5.25) 


i?=[0RjbTA:_1£.]  (5.26) 

For  the  optimization  in  Chapter  VI, 

0m  =  0k  =  0Q  (5.27) 
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VI.  Integrated  Structure/Control  Optimization 


This  Chapter  introduces  the  structure/control  optimization 
procedure.  The  effects  of  Q/R  selection  on  the  control  system 
specifications  are  examined.  In  addition,  two  optimization  approaches 
are  presented.  The  first  approach  i=  based  on  minimizing  mass  and  total 
system  energy.  The  second  approach,  Onada's  method,  develops  a 
relationship  between  control  system  energy  and  control  system  mass.  The 
combined  mass  of  the  structure  and  the  control  system  are  then 
minimized.  Finally,  the  ABE  is  optimized  using  the  constraints  of  the 
existing  control  system. 

Optimization  Algorithm 

The  combined  optimization  of  a  structure  and  its  control  system 
can  be  divided  into  three  distinct  phases  as  shown  in  Figure  6.1.  An 
initial  structure  is  specified  to  start  the  first  phase.  A  finite 
element  program  is  used  to  solve  for  the  eigenvalues  and  eigenvectors  of 
the  structure.  With  the  eigenvalues  and  eigenvectors  determined,  the 
plant  matrix  [A]  can  be  constructed  for  the  state-space  equations.  The 
state-space  equations  are  used  in  the  second  phase  to  solve  the  optimal 
control  problem.  Depending  on  the  output  from  the  second  phase,  the 
third  phase  would  use  an  optimization  procedure  to  determine  a 
"direction"  to  iterate  the  first  phase  input  parameters. 

This  thesis  will  not  use  a  true  optimization  procedure.  Instead, 
a  range  of  structural  parameters  will  be  used  in  a  form  of  "sensitivity 
analysis."  The  structural  parameter  which  is  varied  is  the  cross- 
sectional  area.  To  maintain  a  proportional  cross-section,  the  side 
dimension  "b"  is  defined  in  terms  of  "a"  based  on  the  current  ABE.  "a" 
is  then  varied  in  10%  increments.  Table  6.1  lists  the  cross-sectional 
areas  and  the  moment’s  of  inertia  for  the  17  beam  iterations.  Figure 
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Figure  6.1  Structural/Control  Optimization 
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6.2  shows  the  linear  relationship  between  the  structural  mass  and  the 
cross-sectional  area. 


The  finite  element  procedure  outlined  in  Chapter  IV  was  used  for 
each  of  the  design  iterations.  The  output  obtained  from  MSC/NASTRAN 
included : 

1.  The  natural  frequencies. 

2.  The  mass  normalized  eigenvectors. 

3.  The  [k]  matrix. 

4.  The  initial  conditions  subject  to  a  1000  N  force. 


Table  6.1  Structural  Iteration  Properties 


Number 

Area  (m2) 

(m4) 

I,  <m4) 

J  (m4) 

1 

7 . 900E-05 

6.930E-10 

8 . 426E-10 

2 

1 . 234E-04 

9 . 530E-10 

1 . 692E-09 

2 . 057E-09 

3 

1 . 778E-04 

1 . 976E-09 

3 . 508E-09 

4 . 266E-09 

4 

2 . 419E-04 

3 . 661E-09 

6 . 500E-09 

7 . 903E-09 

5 

3 . 160E-04 

6 . 245E-09 

1 . 109E-08 

1 . 348E-08 

6 

3 . 999E-04 

1 . 000E-08 

1 . 776E-08 

2 . 160E-08 

7 

4 . 938E-04 

1 . 525E-08 

2 . 707E-08 

3 . 292E-08 

8 

5 . 975E-04 

2 . 232E-08 

3 . 964E-08 

4 . 819E-08 

9 

7 . 110E-04 

3. 162E-08 

5 . 614E-08 

6 . 825E-08 

10 

8 . 345E-04 

4 . 355E-08 

7 . 732E-08 

9 . 401E-08 

11 

9 . 678E-04 

5 . 857E-08 

1 . 040E-07 

1 . 264E-07 

12 

1 . 111E-03 

7 . 719E-08 

1 . 370E-07 

1.666E-07 

13 

1 . 264E-03 

9 . 993E-08 

1 . 774E-07 

2 . 157E-07 

14 

1 . 427E-03 

1 . 273E-07 

2 . 261E-07 

2 . 749E-07 

15 

1. 600E-03 

1 . 601E-07 

2 . 842E-07 

3 . 455E-07 

16 

1 . 782E-03 

1 . 987E-07 

3 . 528E-07 

4 . 290E-07 

17 

1 . 975E-03 

2 . 440E-07 

4 . 331E-07 

5 . 267E-07 
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The  control  optimization  phase  was  solved  using  MATLAB.  A  MATLAB 
.M  file  (Appendix  D)  was  written  to  apply  Linear  Quadratic  Regulator 
theory  to  the  finite  element  output  data.  The  MATLAB  .M  file  produced 
data  files  containing  the  following  informations 

1.  The  optimum  gain  matrix  K. 

2.  The  Riccati  solution  matrix  S. 

3.  The  cost  function  J. 

4.  The  initial  control  force  for  each  actuator. 

5.  The  control  energy  function  C. 

6.  The  closed-loop  damping  ratios. 

7.  The  closed-loop  frequencies. 

Each  iteration  was  performed  using  Q  and  R  values  as  defined  is  Chapter 
V  as  well  as  constant  Q  and  R  values.  In  addition,  a  weighing  value  was 
applied  to  the  Q  and  R  values  to  determine  the  effect  weighing  has  on 
the  optimization.  The  basic  nomenclature  scheme  used  is 
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Q  ##  =  Variable  Q  as  defined  in  Equation  5.25 

R  ##  =  Variable  R  as  defined  in  Equation  5.26 

QI  ##  =  Constant  Q  as  defined  in  Equation  E.l 

RI  ##  =  Constant  R  as  defined  in  Equation  E.2 

The  first  #  in  the  scheme  represents  the  weight  applied  to  0Q.  The 
second  #  represents  the  weight  applied  to  0R.  The  weighing  #'s  used 
were  combinations  of  1.0  and  0.3.  Table  E.l  in  Appendix  E  provides  the 
weighing  used  in  the  optimization. 

The  final  phase  of  the  combined  optimization  procedure  uses  data 
from  the  second  phase  to  select  a  search  direction  in  which  to  vary  the 
input  structural  parameters.  As  mentioned  earlier,  a  search  direction 
is  not  used.  Instead,  a  range  of  structural  parameters  are  used  to 
chart  the  effects  of  the  various  optimization  methods.  The  data  needed 
from  the  second  phase  depends  on  the  optimization  method  chosen.  In  the 
following  sections,  the  cost  function  J  and  the  control  energy  function 
C  are  used.  Other  methods  contained  in  the  literature  have  demonstrated 
the  use  of  the  initial  control  force,  the  closed-loop  damping  ratios, 
and  the  closed-loop  eigenvalues  to  find  a  search  direction  to  iterate 
the  structural  parameters. 

Relationship  of  0  and  R  Selection  To  Control  Response 

To  show  the  effects  of  Q  and  R  selection  on  the  optimization 
process,  viscous  damping  from  the  structure  and  actuators  will  be 
neglected.  With  viscous  damping  neglected,  only  damping  from  the 
velocity  feedback  control  law  will  effect  the  system.  Figures  6.3  and 
6.4  show  the  simulated  undamped  response  of  beam  iteration  #1  and  #17 
to  the  applied  impulse  force.  The  measurements  are  taken  from  positions 
11a,  lib,  and  5  on  the  beam. 

From  the  velocity  feedback  control  law,  damping  is  added  to  the 
closed  loop  response.  The  amount  of  damping  depends  on  the  gain  K  which 
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is  selected  by  LQR  theory.  The  closed-loop  plant  matrix  from  equation 
5.14  is  repeated  here. 

Acl  =  [A  -  B*K)  (6.1) 

The  effect  of  selecting  Q  and  R  by  Equations  5.24  and  5.25  on  control 
theory  is  to  keep  the  closed-loop  damping  ratio  and  the  percent 
overshoot  constant  throughout  the  structure/control  iteration.  Figures 
6.5  and  6.6  show  the  closed  loop  damping  ratios  and  the  percent 
overshoot  for  Q  and  R  selected  by  Equations  5.24  and  5.25.  Figures  6.7 
and  6.8  show  the  damping  ratios  and  the  percent  overshoot  when  Q  and  R 
are  held  constant.  From  these  figures,  it  is  easy  to  see  how  the 
weighing  parameters  0  on  Q  and  R  can  be  selected  to  achieve  a  desired 
control  response.  If  we  plow  che  damped  closed-loop  frequency  versus  a 


Figure  6.3  Ream  #1  Undamped  Response 
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Figure  6.4  Beam  #17  Undamped  Response 

=  Cen  ,  we  see  in  Figure  6.9  that  when  the  Q  and  R  values  vary,  a 
negatively  sloped  line  is  produced  in  the  s-plane  which  is  indicative  of 
maintaining  a  constant  overshoot  requirement  (21:61).  A  constant 
settling  time  requirement,  on  the  other  hand,  would  produce  a  vertical 
line  in  the  s-plane.  With  Q  and  R  held  constant.  Figure  6.10  shows  a 
slope  slightly  right  of  a  vertical  line. 

The  settling  time  can  be  found  from  (21:61) 


for  a  settling-time  specification  of  1%.  The  settling  time  is  plotted 
in  Figure  6.11  for  Q  and  R  which  vary  and  Figure  6.12  for  Q  and  R  held 
constant.  The  settling  time  starts  out  high  but  then  decreases  for  the 
Q  and  R  which  vary.  The  settling  time  gradually  increases  for  the 
constant  Q  and  R. 
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Figures  6.13  -  6.16  show  the  initial  control  force  required  from 
each  actuator  for  each  iteration.  In  order  to  maintain  a  constant 
percent  overshoot,  more  force  is  needed  as  the  structure  gets  heavier. 
This  is  shown  for  actuator  A  in  Figure  6.13  and  for  actuator  C  in  Figure 
6.14.  When  a  "constant”  settling  time  is  desired,  more  control  force  is 
needed  as  the  structure  gets  lighter.  This  is  shown  in  Figure  6.15  for 
actuator  A  and  in  Figure  6.16  for  actuator  C. 


Optimization  Using  Cost  Function  J 

The  optimization  of  space  structures  can  be  formulated  (8:63)  as 

minimize  m(v) 

such  that  gj(v)  >0,  j  =  1,  .  .  .  ,  ng 

where  [v]  is  a  vector  of  structural  sizes,  and  m  is  the  structural 
objective  function,  typically  mass.  The  function  g^(v)  represents 
constraints  on  the  structural  design.  To  integrate  structural 
optimization  with  the  control  system  design,  the  mass  objective  function 
can  be  combined  with  the  control  quadratic  performance  index  J,  so  that 
the  optimization  becomes  (3:928) 

minimize 


m*  =  q^w(v)  +  q2J{K,v) 


(6.3) 


such  that  g; ( v)  >  0,  j  =  1,  . . . , ng 

h ,  ( K , v )  >  0 ,  i  =  1,  •  •  •  ,nh 

where  h;  are  constraints  on  the  closed-loop  control  system.  For 
numerical  optimization,  the  gradient  of  m*  is  required.  Miller  and  Shim 
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(22:292)  present  a  numerical  optimization  solution  technique  to  solve 
equation  3. 

Values  of  J  versus  cross-sectional  area  are  plotted  in  Figures 
6.17  -  6.20  for  the  different  weighing  parameters  0  on  Q  and  R.  Figure 
6.17  shows  the  J  values  plotted  for  Q  and  R  selected  by  equations  5.24 
and  5.25.  These  J  values  were  computed  without  viscous  damping  (£  =  0). 
To  show  the  effects  of  damping  on  the  J  values,  Figure  6.18  uses  the 
same  Q  and  R  values  but  includes  viscous  damping.  For  the  lighter 
structures,  the  total  minimum  energy  J  increases  significantly.  Figure 
6.19  is  a  plot  of  J  for  the  undamped  structure  when  Q  and  R  are  held 
constant.  Figure  6.20  is  the  same  plot  except  the  effects  of  damping 
are  included.  When  Q  and  R  are  held  constant,  damping  inherent  to  the 
structure  does  not  have  much  effect  on  the  minimum  system  energy  as 
shown  in  Figure  6.19  and  Figure  6.20.  The  reason  for  this  can  be  seen 
in  Figures  6.15  and  6.16.  For  fixed  Q  and  R,  the  control  system  expends 
more  energy  to  quickly  dampen  out  vibrations  for  the  lighter  structures. 
Even  though  more  control  energy  is  expended,  the  total  system  energy  is 
less  because  the  control  system  does  not  allow  the  structure  to  vibrate 
long  enough  to  dominate  the  total  energy. 

From  equation  6.3,  the  optimal  value  for  m*  is  clearly  dependent 
on  the  choice  of  q1  and  q2.  How  do  we  select  q1  and  q2  when  J(K,v) 
represents  total  system  energy  and  m(v)  represents  structural  mass? 
Clearly  the  units  of  energy  and  mass  are  incompatible,  but  the 
literature  avoids  this  question.  Presumably  in  the  design  process,  q1 
and  q2  would  be  design  values  which  are  adjustable.  q2,  for  example, 
could  be  a  measure  of  mass/energy  for  a  given  control  system.  On  the 
other  hand,  it  may  not  matter.  Since  this  is  multiobjective 
optimization,  q1  and  q2  are  weighing  parameters  which  may  simply  be 
selected  by  decision  makers  to  achieve  some  overall  goal. 

Figure  6.21  and  Figure  6.22  show  the  optimization  of  m(v)  and 
J(K,v)  for  (q,  =  1,  q2  100),  and  (q1  =  1,  q2  =  385),  respectively.  For 


6-9 


q1  =  q2  =  1.0,  the  mass  dominates  the  optimization  and  the  curve  slopes 
up  to  the  right  with  no  minimum.  For  q1  =  1.0  and  q2  >  1000,  the 
control  system  dominates  and  the  curve  slopes  down  the  right  with  no 
minimum.  Somewhere  in  between  the  two  extremes  is  an  optimal  value  for 
q^  and  q2  of  whose  selections  are  another  nested  optimization  problem. 
q2  =  385  was  selected  to  show  a  trend.  As  q2  is  increased,  the  minimum 
m*  moves  to  the  right.  The  same  observation  can  be  made  on  Figures  6.23 
and  6.24  except  that  a  higher  value  of  q2  is  required  for  the  system 
energy  to  influence  the  curve  minimum. 


Optimization  Using  Onada's  Formulation 

Haftka  (8:64)  points  out  two  shortcomings  of  the  above  procedure. 
First,  the  effect  of  the  control  system  on  the  mass  of  the  structure  and 
the  structural  const.aints  is  neglected.  Second,  the  Q  and  R  matrices 
are  somewhat  arbitrary  and  are  tuned  by  the  control  engineer  to  achieve 
the  desired  dynamic  response.  It  is  therefore  unreasonable  to  select 
them  before  the  optimization  procedure  starts  and  then  leave  them 
unchanged.  Problems  with  the  meaning  of  q2  in  the  equation  for 
minimizing  the  mass  and  total  system  energy  were  discussed  above. 

Onada's  formulation  provides  a  physical  meaning  in  a  format  which  can  be 
modeled.  By  selecting  Q  and  R  from  equations  5.25  and  5.26  which  allow 
them  to  vary,  the  second  shortcoming  that  Haftka  points  out  with  the 
above  procedure  is  satisfied. 

Onada’s  formulation  assumes  (8:64)  that  the  mass  of  the  control 
system  is  related  to  the  control  effort  as 

mc  =  a  [  CE  ]  P  (6.4) 

where  a  and  P  are  constants  and  CE  is  the  quadratic  measure  of  the 
control  effort 
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(6.5) 


CE  =  u  0TRu0 

The  optimization  problem  is  to  minimize  the  total  mass  of  the  system 
subject  to  a  constraint  on  the  magnitude  of  the  response.  The  total 
mass  of  the  system  is  the  sum  of  the  control  mass  and  the  structural 
mass.  The  response  of  the  system  is 

r  =  x0tQx0  (6-6) 

The  optimization  problem  is  formulated  as: 

Find  [v]  and  [K]  to  minimize 

m=ms(v)+mc(v,K)  (6.7) 

subject  to  the  constraint  r[v,K]  <  a  where  a  is  the  response  magnitude 
allowable. 

For  a  beam  structure,  Onada  and  Haftka  showed  (23:1136)  that  as 
the  mass  of  the  beam  was  decreased,  the  mass  of  the  control  system 
became  very  large  in  order  to  control  the  vibration  modes  in  the  lighter 
structure.  This  was  also  demonstrated  for  the  ABE  optimization  when  the 
values  of  Q  and  R  were  held  constant. 

Figure  6.25  is  a  plot  of  the  C  values  for  the  undamped  structure 
when  Q  and  R  vary.  The  energy  oi  the  control  system  remains  fairly 
constant.  In  Figure  6.17  then,  the  sharp  increase  in  J  for  the  lighter 
structures  can  be  contributed  to  the  increase  of  structural  kinetic  and 
strain  energy.  Figure  6.26  shows  the  effect  of  damping  on  the  C  values. 
When  Q  and  R  are  held  constant,  the  control  energy  increases  as  the 
structural  size  decreases.  This  is  shown  in  Figures  6.27  and  6.28. 

Once  again,  structural  viscous  damping  does  not  effect  C  significantly. 

Figures  6.29  -  6.32  show  the  plot  of  equation  6.7  where  mc(v,K)  is 
determined  by  equation  6.4  for  a  =  50  and  (5  =  1.  Figures  6.25  and  6.26 
are  for  Q  and  R  which  vary.  Since  C  remains  fairly  constant  when  Q  and 
R  vary,  the  optimal  choice  for  m  is  the  smallest  structural  size  which 
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satisfies  any  given  structural  constraints.  When  Q  and  R  remain 
constant,  however,  C  becomes  large  as  the  structural  size  of  the  beam  is 
decreased.  The  optimization  formulation  of  equation  6.7  produces  a 
minimum  on  Figures  6.31  and  6.32  for  a  =  50.  The  total  mass  of  the 
system  increases  significantly  for  the  smaller  structural  sizes  as  shown 
by  Onada  and  Haftka.  This  is  reasonable  since  a  larger  control  would  be 
required  to  provide  the  higher  energy  levels. 

Comparison  of  Optimization  Methods  and  O/R  Selection 

The  direct  relationship  shown  in  Figures  6.5  -  6.14  to  control 
specification  parameters  gives  an  indication  of  how  to  select  Q  and  R 
and  their  weighing  parameters.  To  maintain  constant  overshoot  and 
closed-loop  damping  requirements,  Q  and  R  are  chosen  by  Equation  5.24 
and  5.25.  The  weights  on  Q  and  R  are  varied  to  select  the  control 
response  desired.  For  a  near  constant  settling  time,  Q  and  R  can  be 
fixed  with  minor  adjustments  to  their  weighing  parameters. 

The  effect  on  the  control  system  is  shown  in  Figures  6.33  to  6.37. 
Comparing  beam  #1  in  Figures  6.33  to  6.34,  more  initial  actuator  force 
is  required  when  Q  and  R  are  held  constant  to  achieve  a  constant 
settling  time.  When  Q  and  R  vary,  the  initial  actuator  force  is  smaller 
but  the  time  interval  over  which  it  is  applied  is  longer.  For  actuator 
A,  the  square-root  of  the  squared  forces  for  both  Q/R  selections  are 
plotted  in  Figure  6.35  to  compare  the  force  magnitudes  with  their 
applied  time  intervals.  As  the  structure  gets  larger,  more  control 
force  is  required  to  maintain  a  constant  overshoot  and  closed-loop 
damping  ratio.  Figure  6.36  shows  this  increase  in  initial  force  when  Q 
and  R  are  varied.  For  the  larger  structure,  less  initial  control  force 
is  needed  to  achieve  a  constant  settling  time  as  shown  in  Figure  6.37. 

It  is  interesting  to  examine  the  structural  response  to  these 
control  inputs.  Figures  6.38  and  6.40  compare  beam  #1  and  beam  #17  when 
Q  and  R  are  varied.  In  Figure  6.40,  the  structural  response  is  damped 
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out  very  quickly  to  maintain  the  overshoot  requirements.  The  response 
in  Figure  6.38  is  allows  to  vibrate  for  some  time  (3+  sec)  because  the 
overshoot  requirement  is  not  violated.  Figures  6.39  and  6.41  compare 
beam  #1  and  beam  #17  when  Q  and  R  are  held  constant.  In  both  cases,  the 
vibration  is  damped  out  in  less  than  1  sec.  Less  control  input  is 
required  to  achieve  this  response,  however,  because  of  the  increased 
damping  present  in  the  larger  structure. 

It  is  also  of  interest  to  compare  the  J  and  C  values  for  the 
different  Q  and  R  selections.  Figures  6.42  and  6.43  are  plots  of  J 
versus  area  and  C  versus  area,  respectively,  for  both  choices  of  Q  and 
R.  The  total  system  energy  starts  out  lower  when  Q  and  R  are  fixed  than 
when  Q  and  R  are  varied.  At  some  point  on  the  curve,  the  total  system 
energy  becomes  less  when  Q  and  R  are  varied.  The  same  pattern  emerges 
in  Figure  6.43  when  we  compare  C  values  except  that  the  trend  is 
reversed  for  the  Q  and  R  selection.  Notice  that  the  cross-over  point 

for  J  and  C  values  are  not  the  same.  At  the  point  of  cross-over  on  the 
J  curve,  our  expected  value  for  C  is  less  for  the  case  when  Q  and  R  are 
constant . 

Selection  of  Structure  By  Sensitivity  Analysis 

A  true  optimization  was  not  performed  on  the  ABE,  but  we  can  use 
the  trends  presented  in  this  Chapter  to  examine  the  current  ABE 
structure.  Since  the  values  of  a  and  P  are  not  known  for  the  ABE,  the 
minimization  of  the  total  system  energy  will  be  used. 

From  Chapter  III  and  Appendix  B,  the  limits  on  the  proof  mass 
actuators  and  the  structural  dynamics  shaker  can  be  found.  These  limits 
are 

Actuators:  8.9  N  *  80%  =  7.12  N  each 

Shaker:  45  N 
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In  addition,  the  actuators  are  limited  to  the  bandwidth  0-100  Hz. 

Figure  6.44  shows  the  initial  control  force  required  for  actuator  A  (and 
B)  for  the  different  parameters  8  when  Q  and  R  are  varied.  Figure  6.45 
shows  the  initial  control  force  required  for  actuator  C  under  the  same 
conditions.  A  constraint  line  based  on  the  above  limits  is  drawn  on  the 
diagrams . 

The  plot  of  equation  6.3  with  q1  =  1.0  ana  q2  =  385.0  is  shown  in 
Figure  6.46.  For  the  top  three  curves,  there  are  two  minimums  on  the 
graph.  Which  minimum  to  select  depends  on  the  limits  of  our  control 
system.  From  Figures  6.44  and  6.45,  we  see  that  the  Q0/R0  and  the  Q0/R3 
curve  require  initial  control  forces  beyond  the  limits  of  our  actuators. 
The  required  initial  control  forces  for  Q3/R0  are  within  our  limits. 
Q3/R0  is  selected  as  the  minimum.  For  the  selection  of  q1  =  1.0  and 
q2  =  385,  the  minimum  for  happens  to  coincide  with  the  existing  ABE 
structure.  A  plot  of  the  required  controller  response  is  shown  in 
Figure  6.47  and  a  plot  of  the  structural  response  is  shown  in  Figure 
6.48.  Of  course  other  fl  and  q,-  parameters  for  the  ABE  will  also  produce 
minimums  within  the  limits  of  our  controllers.  These  parameters  can  be 
adjusted  to  determine  a  control  law  which  produces  the  desired  system 
response  with  the  limits  of  the  control  system.  The  control  law  for 
Q3/R0  is  presented  in  Appendix  F. 


6-14 


0.35 


QIO.RIO' 

QI3.RI0' 


6-17 


Figure  6.7.  Damping  Ratios  -  Q/R  Constant  (Undamped  Structure) 
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Figure  6.8.  Percent  Overshoot  -  Q/R  Constant  (Undamped  Structure) 
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Figure  6.9.  S-Plane  -  Q/R  Vary  (Undamped  Structure) 
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Figure  6.10.  S-Plane  -  Q/R  Constant  (Undamped  structure) 
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.12  Settling-Time  -  Q/R  Constant  (Undamped  Structure) 
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Figure  6.13  Actuator  A  Initial  Force  -  Q/R  Vary  (Undamped  Structure) 
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Figure  6.14  Actuator  C  Initial  Force  -  Q/R  Vary  (Damped  Structure) 
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15  Actuator  A  Initial  Force  -  Q/R  Constant  (Undamped  Structure) 
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.16  Actuator  C  Initial  Force  -  Q/R  Constant  (Damped  Structure) 
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Figure  6.17  J  Values  -  Q/R  Vary  (Undamped  Structure) 
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Figure  6.18  J  Values  -  Q/R  Vary  (Damped  Structure) 
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Figure  6.19  J  Values  -  Q/R  Constant  (Undamped  Structure) 
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Figure  6.20  J  Values  -  Q/R  Constant  (Damped  Structure) 
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Figure  6.21  Equation  6.3  -  Q/R  Vary  (q1  =  1.0;q,  =  100.0) 


.22  Equation  6.3  -  Q/R  Vary  (q1  =  1.0;q2  =  385.0) 
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Figure  6.24  Equation  6.3  -  Q/R  Constant  (q1  =  1.0;q2  =  385.0) 


0  30 


6-35 


Figure  6.25  C  Values  -  Q/R  Vary  (Undamped  Structure) 
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Figure  6.27  C  Values  -  Q/R  Constant  (Undamped  Structure) 
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Figure  6.28  C  Values  -  Q/R  Constant  (Damped  Structure) 
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Figure  6.29  Equation  6.7  -  Q/R  Vary  (Undamped  Structure, 


Figure  6.30  Equation  6.7  -  Q/R  Vary  (Damped  Structure, 
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Figure  6.31  Equation  6.7  -  Q/R  Constant  (Undamped  Structure, 


00  93 


6-42 


Figure  6.32  Equation  6.7  -  Q/R  Constant  (Damped  Structure, 
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Figure  6.36  Beam  #17  Actuator  Responses  -  Q/R  Vary  (Undamped  Structure  0Q=8R=1) 


6-47 


Time  ( s e c  ) 


r; 

i 


6-48 


Figure  6.38  Beam  #1  Structural  Response  -  Q/R  Vary  (Undamped  Structure 
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Figure  6.39  Beam  #1  Structural  Respbnse  -  Q/R  Constant  (Undamped  Structure 
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Figure  6.40  Beam  #17  Structural  Response  Q/R  Vary  (Undamped  Structure 
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Figure  6.42  J  Value  Comparison 
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Figure  6.43  C  Value  Comparison 
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Figure  6.44  Initial  Force  For  Actuators  A  &  B  -  Damped  Structure 
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Figure  6.45  Initial  Force  For  Actuator  C  -  Damped  Structure 


QO.RO 


6-56 


Figure  6.46  Optimization  of  Damped  Structure  Using 
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VII.  Conclusions  and  Recommendations 


Conclusions 

Our  goal  is  to  place  large  structures  in  space.  Both  NASA  and  DOD 
are  working  on  systems  to  be  launched  by  the  year  2000.  If  we  are  to 
realize  this  goal,  we  must  explore  every  avenue  to  limit  the  economic 
costs.  The  cost  of  launching  large  structures  into  space  will  dominate 
the  overall  system  cost.  Since  launch  costs  are  driven  by  weight, 
reducing  the  weight  of  our  structures  will  help  keep  them  within  our 
economic  limits. 

The  purpose  of  this  work  was  to  model  an  existing  structure  and 
then  to  examine  the  structure/control  optimization  process.  Since 
weight  reduction  was  desired,  two  optimization  processes  were  selected 
which  minimize  mass.  The  first  method  optimized  the  combined  mass  of 
the  structure  and  the  total  system  energy  of  the  structure  and  the 
control  system.  The  second  method  related  the  mass  of  the  control 
system  to  the  required  control  energy,  then  optimized  the  mass  of  the 
control  system  and  the  mass  of  the  structure.  Both  methods  proved 
successful,  but  each  method  resulted  in  different  optimized  structures. 

The  reason  for  the  difference  can  be  linked  to  the  selection  of  Q 
and  R.  When  Q  and  R  are  allowed  to  vary,  the  effect  on  control  theory 
is  to  maintain  a  constant  percent  overshoot.  For  the  lighter 
structures,  a  smaller  control  system  can  prevent  the  system  from 
exceeding  the  overshoot  specification.  However,  the  smaller  control 
system  will  not  dampen  the  vibration  as  quickly.  Since  the  first  method 
uses  the  total  system  energy,  a  higher  J  value  will  result  for  the 
smaller  structures.  The  goal  is  to  then  lower  this  system  energy  by 
increasing  the  structural  mass  and  the  control  system  inputs  to  a  point 
where  system  mass  and  system  energy  are  optimized.  Onada's  formulation 
has  the  opposite  effect.  Since  the  control  effort  required  is  lower 
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when  constant  overshoot  requirements  are  needed,  the  optimization 
implies  that  the  smallest  structure  which  meets  the  load  specifications 
is  the  optimized  structure.  When  settling  time  requirements  drive  the 
control  design,  more  control  energy  is  needed  to  quickly  dampen  out  the 
lighter  structures.  This  results  in  a  higher  C  value  for  Onada's 
formulation  but  since  the  vibrations  have  been  quickly  attenuated,  the 
overall  J  value  is  small.  Onada’s  optimization  results  in  very  heavy 
control  systems  then  for  the  smaller  structures. 

Contributions 

Which  optimization  method  do  you  use  then  and  how  do  you  select  Q 
and  R  with  their  appropriate  weighing  parameters?  The  answer  to  this 
question  is  highly  dependent  on  the  design  goals  you  are  trying  to 
achieve.  Trends  were  shown  which  should  help  other  researchers  select 
Q/R  and  an  appropriate  optimization  procedure  for  various  control 
specifications.  For  the  ABE,  this  work  has  provided  a  baseline  for 
additional  research  on  the  structure/control  optimization  problem. 

Recommendations  For  Future  Work 

1.  Additional  work  in  needed  to  characterize  the  actuator 
dynamics.  The  finite  element  model  created  for  this  work  included  the 
viscous  damping  of  the  actuators.  The  viscous  damping  assumption  was 
verified  (7:19)  by  modal  testing.  Stiffness  effects  from  the  electronic 
centering  spring,  however,  were  not  modeled  and  this  contributed  to  the 
high  error  rate  in  Table  4.6. 

2.  One  of  the  biggest  challenges  (4:472)  facing  LSS  control 
engineers  is  the  development  of  inertial  actuators.  The  ABE  provides  an 
excellent  test  bed  on  which  to  design  new  actuators  and  to  compare  the 
performance  of  different  actuators.  Cristler  (5)  goes  to  some  length  in 
documenting  the  problems  associated  with  developing  the  existing  proof 
mass  ac’-uators.  Theoretical  work  charting  the  stability  regions  of  the 
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structure/control  system  would  help  future  engineers  design  new 
actuatorj.  This  work  should  include  the  effects  of  damping,  stiffness, 
and  control 

3.  To  examine  the  control/structure  optimization  process,  this 
won  assumed  that  all  modes  were  observable.  Clearly  this  is  not  the 
case.  The  initial  conditions  suggest  that  higher  order  modes  are  being 
excited  to  a  greater  degree  in  the  smaller  structures  than  in  the  larger 
structures.  A  comparison  between  the  results  of  this  work  and  one  which 
in;l”dos  a  state  estimator  woul  ’  be  interesting  to  determine  the  effect 
higher  order  modes  have  on  th*  selection  of  the  optimal  structure- 
control  sjsc^n’. 

4.  A  "true"  optimum  was  not  determined  in  this  work.  An 
opportunity  exists  to  program  the  methods  presented  in  Chapter  VI 
subject  to  the  constraints  :n  the  ABE  control  system.  In  addition,  the 
selection  of  q1  and  q2  can  be  nested  in  the  optimization  process.  If 
physical  meaning  can  be  given  to  q1  and  q2,  their  inclusion  in  the 
optimization  loop  will  be  necessary  to  achieve  a  true  optimum  design. 

5.  Damping  inherent  in  the  both  the  structure  and  the  control 
system  was  modeled  by  proportional  damping.  For  the  lighter  structures, 
inertia  played  a  greater  role  in  the  damping  factor.  Stiffness,  on  the 
other  hand,  was  predominate  for  the  larger  structures.  As  a  parallel, 
different  parameters  were  used  to  weight  Q  and  R  for  comparison.  When  Q 
and  R  were  selected  by  equations  5.25  and  S."^,  it  was  assumed  that  0m  = 
0k.  What  effect  does  varying  these  parameters  have  on  the  structure- 
control  optimization? 

6.  Finally,  Jacques  pointed  out  in  his  conc'usion  (7:87-88)  that 
the  ABE  was  an  excellent  test  bed  for  further  controls  experiments  and 
research.  Additional  work  can  be  done  to  compare  the  various  control 
laws  and  estimator  designs.  Other  optimization  schemes  can  be  tested  on 
the  ABE  before  their  use  on  more  complicated  structures.  In  short,  new 
ideas  should  be  tested  in  a  known  environment  so  their  true  merits  can 
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be  determined.  With  the  base  line  experimental  and  theoretical  research 
accomplished,  the  ABE  is  a  known  environment. 
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Appendix  A:  Calculation  of  ABE  Physical  Properties 


Beam  Physical  Properties 


I  =ab 

V  ~TT 


T  _  a3jb 

Z~TT 


Jp=Iy+Iz 


(A.l) 


(A. 2) 

(A. 3) 


A-l 


In  general,  if  the  cross-section  is  not  circular,  there  will  be  some 
warping  of  the  cross-sectional  plane.  For  circular  cross  sections 
(9:200) , 


For  noncircular  cross  sections,  however. 


J  can  be  obtained  from  texts  on  advanced  strength  of  materials. 
(24:290),  the  value  of  J  for  a  rectangular  beam  is 
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Disk  Physical  Properties 


(A. 4) 


(A. 5) 
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Figure  A. 2  Disk  Dimensions 
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DEFINITIONS  (Continued) 
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rat-ng  can  oe  ca>c-.atec  'rQr  P«-2  A  me  resistance  usee  -n  m*s  cj'Cuia*<on  it  •<  miiimuw  temoeratu'«  not  me  *oo m  temoe'a’jre  .i  ^e 

g>v«"  on  r r>t  a»u  snee‘ 

Tn«  'essences  o<  do  m  me  con  anc  m«  secondary  increase  at  0  **«  oe'  a  eg'**  Cei$>us  ana  me  tytxca1  temperature  r.*e  oi  me  seconaa'y  •* 
aooui  ;  3  mat  o'  me  co*'  ams  current  rating  s  caicuiatea  "om  tr»«*  iniofmat-or  O'  measu'eo  eipe'imentatiy 

AMS  current  rating  muitipnefl  Dy  »Orct  constant  gives  AMS  'ore*  rating  Tpms  it  a  ute'ui  numo*t  of  compering  caoaO>'>t'et  0*  motors  ana  It  Sim.i»f 
to  ratec  toraue  o'  a  rotary  moto' 


ngone  The  ot*e»  t»o  orocesses  oeoena  on 
whoee  iu*  ooootet  armature  *u,a  Decays 
t  given  on  me  oa*a  sneet  it  a  conservative 


Figure  B.l  Linear  Motor  Specifications 
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HPD  SERIES— HERMETICALLY  SEALED  (PIN  TERMINATION) 


•  HERMETICALLY  SEALED  BY  TIG  AND 
EB  WELDING 

•  IMPERVIOUS  TO  HOSTILE  ENVIRONMENTS 

•  THROUGH-BORE  CONSTRUCTION 

HPD  Scntt  unitsere  similar  to  the  DC-D  end  HCD  Senes. 
Tungsten  men  gas  (TIG)  and  electron  beam  (EBI  weld¬ 
ing  provide  hermetic  sealing  that  is  free  from  oxidation- 


GENERAL  SPECIFICATIONS 

Inout . *15  V  DC  Inommall.  ±20  mA 

Operating  Temper¬ 
ature  Range . 0*F  to  *160*F  l-18*Cto  *70*CI 

Survival  Temper¬ 
ature  Range . -85*F  to  *200*  F 

1-55‘C  to  -95*0 

Null  Voltage  . 0  V  DC 

Ripple  . Less  than.  25  mV  rms 

Linearity . *  0.25%  full  range 

Stability . 0.125%  full  scale 


producing  faults  that  may  cause  leakage .  For  this  reason. 
HPD  Sents  LVDT's  are  impervious  to  dirt,  water,  steam 
spray,  and  most  corrosives.  They  have  been  Qualified  et 
pressures  up  to  1000  psig  (70  bars)  and  are  suitable  tor 
numerous  high-pressure  applications.  HPD  units  employ 
a  glass-sealed,  pm-terminal  header  that  allows  the  core 
end  core  rod  to  pass  through  the  unit.  HPD  units  have 
oouDle  magnetic  shielding  that  makes  them  insensitive 
to  external  magnetic  influences. 


Temperature 
Coefficient  of 
Scale  Factor  .  . 
Shock  Survival  .  .  . 
Vibration  Toleranoe 
Coil  Form  Material 
Housing  Material  . 
Electrical 

Termination  .  . 
Output  impedance 


0.04%/*F  (0.08  VC  I 

250  g  for  11  milliseconds 

10;  up  to  2  kHz 

High  Density,  glass-filled  polymer 

AISl  400  senes  stainless  steel 

6-pm  terminal  header 
Less  than  1  Onm 


PERFORMANCE  SPECIFICATIONS  AND  DIMENSIONS 
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Figure  B.3  LVDT  Specifications 
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APS  113  with 
0112  Reaction  Mae 


APS  111  Atr  ••■tins 


OPTIONAL  CONFIGURATIONS 

APS  111- AS  Air  Seating  Modal  Atr 

lubricated  Bushings  replace  ina  linear 
Mil  Bushings  used  m  ina  Oik 
ELECTRO-SEIS  armatura 
guidance  sysrem  in  addition  an  ax 
distribution  system.  tia  down  and 
leveling  Base  aia  provided 

The  naar  zero  iriciion  o<  me  an 
Bushings  i  an  essential  feature  'o' 
maaiu'ing  resonance  decay  ales  in 
very  hgntiy  damped  structures 

Tha  An  Baa'mg  coniigu'alion  estends 
tha  aooi'caiion  of  tha  Banc  APS  til  10 
include  tha  calibration  and  actuation 
o'  acceiaiomaiau  and  othat  motion 
transducers  m  ina  teurmc  traguancy 
'anga 


APS  113-LZ  Lota  Impedance  Cod 
All  features  of  the  Banc 
ELECTRO-SEIS  Shaker  are 

retained  The  drive  coil  it  wound  m  a 
manner  which  allow  lanes  or  parallel 
connection,  offering  tne  user  the 
choice  of  standard  or  low  impedance 
Thu  cpnon  is  required  it  the  shaker  is 
to  Be  used  with  the  APS  12* 

SuAL  bSCCC  Power  Amplifier  for 
axended  freouency  ’ange  or  random 
noise  aicitaiion 


APS  111-HP  High  Perce  Cod 

aii  features  of  me  oasic 
ELECTRO-EE  IE  Shaker  are 
'etamed  as  m  the  APS  1 13-LZ  The 
drive  con  is  provided  to  match  me 
APS  <2*  SUAL-CSSC  Power 
Amplifier  tor  *0%  increase  m  lorce 
with  a  SON  duty  cycle  1 1. 2  hr  cycle) 


APS  1 13- LA  Lightweight  Armature 
The  Body  ot  the  ELECTRO-SEIS 
Shaker  is  retained  But  tne  armature 
and  guidance  system  are  'ecaced  with 
elements  ottering  substantial  weight 
reduction  The  onve  cou  is  iighiened 
—  with  corresponding  reduction  m 
masurium  force  —  ana  the  armature 
guidance  sytrem  elements  are  reduced 
m  sue  and  wdgnt  Tn,s  results  m  a 
corresponding  -eduction  m  cross  art 
stiff  ness  and  >oad  carrying  aoiiity  Tne 
long  stroke  capability  •*  'etamed  and 
the  treouency  range  tor  masurium 
torce  output  Ik  extended  to  1000  Hz 

The  Ligntweignr  Armature  -a  a 
desirame  feature  w "»«  using  me 
Shaker  lor  ((Citing  structures  hev.ng 
•ow  mooai  mass 
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Shaker  Specifications 


ACOUSTIC  POWER  SYSTEMS,  INC 

Systems  lor  Generating  Controlled  Vibration 


SPECIFICATIONS 

Mode*  114 

Model  124 

Modal  124- EP 

Average  Output  into  shaker  reactive  load 

125  V-A  rms 

250  V-A  rms 

250  V-A  rms 

pea«  Output,  into  snaker  reactive  load 

250  V-A  rms 

500  V-A  rms 

750  V-A  peak 

Current  Output  peak  iranoom  no>sei 

60  A  peak 

12  0  A  peak 

IS  0  A  peak 

Current  Output  continuous 

40  A  rms 

SO  A  rms 

80  Amts 

Frequency  Range 

0-2000  Ha 

0-2000  hi 

0-2000  Hi 

input  Sk)nai  Voltage 

2  V  peak 

2  V  peak 

2  V  peak 

input  Impedance 

too  K  ohm 

100  K  ohm 

lOOKonm 

Noise  —  reterred  to  man  output 

-90  d6 

-90  dB 

-90  dB 

Current  Moneor  output 

250  mVrA 

125  mV  A 

125  mV  A 

input  Power 

120  V  50  60  Hi.  300  W 
220  2*0  V  optional 

120  V  50  SO  Hi  600  W 
220  240  V  optional 

120  V  50  60  Hi  600  W 
240  V  optional 

Rear  Panel  Connectors 

Power  Output 

WK3-31S  iCannonl 

INK3-31S  iCannont 

VKK3-31S  (Cannon) 

input  Current  Monitor 

BnC  Type  3  aa 

BNC  Type  3  ea 

SNC  Type  3  aa 

AC  Power 

Std  3- Pm  Receptacle 

Sid  3-Pm  Receptacle 

Std  3- Pm  Race  place 

Weigm 

25  ip  ( 11  3  kg| 

45®  120  kg) 

45  «  120  kg| 

Sice  h.in.o 

5  22  «  17  «  9  25  mcnes 

133  ■  432  k  23S  mm 

5  22  i  17.  13  25 mcnes 

133  «  432 . 337  mm 

5  22.17.1325  mcnes 
133 . 432 . 337  mm 

SPECIFICATIONS 

Modal  IIS 

mput  s-gnal  Laval 

t  Vpeaa 

Input  Impedance 

100  K  ohm 

Output  Sgnai  Level 

0  to  to  V  peak 
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Figure  B.5  Acoustic  Power  Systems  Model  114  Power  Amplifier 
Specifications 


Appendix  C:  MSC/NASTRAN  Input  Decks 


Input  Deck  1:  Orainal  ABE  with  impact  load;  no  damping. 


ASSIGN  OUTPUT2=RUN2 .OP2 , STATUS=NEW, UNIT=12 , FORMAT=UNFORMATTED 
TIME  30 
DIAG  64 

SOL  112  $ MODAL  TRANSIENT  RESPONSE 

CEND 

TITLE  =  ABE  FEM  —  UNITS  =  SI 

METHOD  =97 

DLOAD  =  40 

SET  1  =  5  11 

TSTEP  =  98 

S VECTOR  =  ALL 

VELOCITY  =  1 

DISPLACEMENT  =  1 

SPC  =  1 


SUBCASE 

1 

SUBTITLE  = 

CASE  SET 

1 

BEGIN 

GRID 

BULK 

1 

o 

o 

0. 

0. 

0 

GRID 

2 

01 . 78E-01 

0. 

0. 

0 

GRID 

3 

03 . 56E-01 

0. 

0. 

0 

GRID 

4 

05.34E-01 

0. 

0. 

0 

GRID 

5 

07.12E-01 

0. 

0. 

0 

GRID 

6 

0  .889 

0. 

0. 

0 

GRID 

7 

0  1.067 

0. 

0. 

c 

GRID 

8 

0  1.245 

0. 

0. 

0 

GRID 

9 

0  1.422 

0. 

0. 

0 

GRID 

10 

0  1.6 

0. 

0. 

0 

GRID 

11 

0  1.797 

0. 

0. 

0 

CBAR 

+EA 

1 

1 

1  1 

2 

1. 

1. 

+EA 

1 

0 

CBAR 

+EA 

2 

2 

1  2 

3 

1. 

1. 

+EA 

2 

0 

CBAR 

+EA 

3 

3 

1  3 

4 

1. 

1. 

+EA 

3 

0 

CBAR 

+EA 

4 

4 

1  4 

5 

1. 

1. 

C-l 


+EA 

4 

0 

CBAR 

+EA 

5 

5  1 

5 

6 

1.  1. 

+EA 

5 

0 

CBAR 

+EA 

6 

6  1 

6 

7 

1.  1. 

+EA 

6 

0 

CBAR 

+EA 

7 

7  1 

7 

8 

1.  1. 

+EA 

7 

0 

CBAR 

+EA 

8 

8  1 

8 

9 

1.  1. 

+EA 

8 

0 

CBAR 

+EA 

9 

9  1 

9 

10 

1.  1. 

+EA 

9 

0 

CBAR 

+EA 

10 

10  1 

10 

11 

1.  1. 

+EA 

10 

0 

CONM2 , 

11, 

11,0,4.986,0. ,0. ,0. , , 

+EA11 

+EA11 , 

5 . 79E-02 , 0 . , 2 . 895E-2 , 0. , 

0. ,2.895E 

-2 

CONM2 , 

12, 

11,-1, 2. 49E-02, 1.7755 

,-.0572,- 

.1016, 

, +EA12 

+EA12, . 17E-05 , 0 . , .55E-05,0. ,0.  , .40E-05 

CONM2,13f 11 ,-1,3. 49E-02,1. 7678,0. ,-.1016, ,+EA13 

+EA13 , 3 . 81E-05 , 0 . , . 018E-05 , 0 . , 0 . , 3 . 81E-05 

CONM2 ,14,11,-1,8. 48E-02 ,1.7755, . 0572 ,-. 1016, , +EA14 

+EA14 , . 59E-05 , 0 . , 1 . 87E-05 , 0 . , 0 . , 1 . 36E-05 

CONM2, 15, 11,-1, 3. 49E-02, 1.727, 0.0572,-. 1016, ,+EA15 

+EA15 , . 24E-05 , 0. , 1 . 07E-05 , 0. , 0 . , .86E-05 

CONM2 ,16, 11, -1, 9. 47E-02, 1.7132, .062,-. 1016, ,+EA16 

+EA16, 15 . 62E-5 , 0 . , . 43E-05 , 0 . , 0 . , 15 . 62E-5 

CONM2 ,17,11,-1,2. 49E-02 ,1.7755, .0572, .1016, ,+EA17 

+EA17 , . 17E-05 , 0 . , . 55E-05 , 0 . , 0 . , .40E-05 

CONM2 ,18,11,-1,3. 49E-02 , 1 . 7678,0 . , . 1016, , +EA18 

+EA18 , 3 . 81E-05 , 0 . , . 18E-06, 0 . , 0 . , 3 . 81E-05 

CONM2,19, 11, -1,8. 48E-02, 1.7755, -.0572, . 1016, ,+EA19 

+EA19, . 59E-05 , 0 . ,1.87E-05,0. ,0. ,1.36E-05 

CONM2, 20, 11,-1, 3. 49E-02, 1.727, -.0572, .1016, ,+EA20 

+EA20, . 24E-05 , 0 . , 1 . 07E-05 , 0 . , 0 . , .86E-05 

CONM2, 21, 11,-1, 9. 47E-02, 1.7132, -.062, .1016, ,+EA21 

+EA21 , 15 . 62E-5 , 0 . , . 43E-05 , 0 . , 0 . , 15 . 62E-5 

CONM2, 22, 11,-1, 51. 36E-2, 1.7678,0. ,-.1016, ,+EA22 

+EA22 , 40 . 99E-5 , 0 . , 33 . 13E-5 , 0 . ,0. ,40.99E-5 

CONM2 ,23,11,-1,51. 36E-2 ,1.7678,0. , .1016, ,+EA23 

+EA23 , 40 . 99E-5 , 0 . , 33 . 13E-5, 0. ,0. ,40.99E-5 

MAT1  17 . 45E+102 . 83E+10  .317  2766 . 3 . 61E-06 

+MA  1 

+MA  1  1.5E+09  1.5E+09  6.8E+07 

PBAR  1  14.94E-042. 71E-081 . 53E-083 . 30E-08  0. 

+PA  1 

+PA  1  0.  0.  0.  0.  0.  0.  0. 

O.+PB  1 

+PB  16 . 68E-011 . 19E+00  0. 
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SPC 


1 


1  123456 


0. 


SPC  1 

2 

135 

0. 

SPC  1 

3 

135 

0. 

SPC  1 

4 

135 

0. 

SPC  1 

5 

135 

0. 

SPC  1 

6 

135 

0. 

SPC  1 

7 

135 

0. 

SPC  1 

8 

135 

0. 

SPC  1 

9 

135 

0. 

SPC  1 

10 

135 

0. 

SPC  1 

11 

135 

0. 

TLOAD1,40,41,  ,0,42 

DAREA, 41,11,2, 1000. 

TABLED 1 , 42 ,,,,,,, , +E42 
+E42 , 0 . ,0. ,0.001,0. , .002, 

1. , .003,0. 

,  +E43 

+E43, .004,0. ,ENDT 
PARAM  AUTOSPC  YES 

EIGR,97,MGIV,0. ,150. 

,3,3, 

,  ,  ABC 

+BC, MASS 

TSTEP,98, 10, .001,1, , , , ,+E98 
+E98,, 100, .01,1 
P ARAM, POST, -2 
ENDDATA 


Input  Deck  2:  Oroinal  ABE  with  impact  load;  damping  included 


ASSIGN  OUTPUT2=RUN2 .OP2 , STATUS=NEW,UNIT=12 , FORMAT=UNFORMATTED 
TIME  30 
DIAG  64 

SOL  112  $ MODAL  TRANSIENT  RESPONSE 

CEND 

TITLE  =  ABE  FEM  —  UNITS  =  SI 

METHOD  =97 

SDAMPING  =  99 

DLOAD  =  40 

SET  1  =  5  11 

TSTEP  =  98 

SVECTOR  =  ALL 

VELOCITY  =  1 

DISPLACEMENT  =  1 

SPC  =  1 

SUBCASE  1 

SUBTITLE  =  CASE  SET  1 
BEGIN  BULK 

GRID  1  0  0.  0.  0.  0 
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GRID 

2 

01. 

78E-01 

0. 

0. 

0 

GRID 

3 

03. 

56E-01 

0. 

0. 

0 

GRID 

4 

05. 

34E-01 

0. 

0. 

0 

GRID 

5 

07. 

12E-01 

0. 

0. 

0 

GRID 

6 

0 

.889 

0. 

0. 

0 

GRID 

7 

0 

1.067 

o. 

0. 

0 

GRID 

8 

0 

1.245 

0. 

0. 

0 

GRID 

9 

0 

1.422 

0. 

0. 

0 

GRID 

10 

0 

1.6 

0. 

0. 

0 

GRID 

11 

0 

1  -  797 

0. 

0. 

0 

CBAR 

+EA 

1 

1 

1 

1 

2 

1. 

1. 

+EA 

1 

0 

CBAR 

+EA 

2 

2 

1 

2 

3 

1. 

1. 

+EA 

2 

0 

CBAR 

+EA 

3 

3 

1 

3 

4 

1. 

1. 

+EA 

3 

0 

CBAR 

+EA 

4 

4 

1 

4 

5 

1. 

1. 

+EA 

4 

0 

CBAR 

+EA 

5 

5 

1 

5 

6 

1. 

1. 

+EA 

5 

0 

CBAR 

+EA 

6 

6 

1 

6 

7 

1. 

1. 

+EA 

6 

0 

CBAR 

+EA 

7 

7 

1 

7 

8 

1. 

1. 

+EA 

7 

0 

CBAR 

+EA 

8 

8 

1 

8 

9 

1. 

1. 

+EA 

8 

0 

CBAR 

+EA 

9 

9 

1 

9 

10 

1. 

1. 

+EA 

9 

0 

CBAR 

+EA 

10 

10 

1 

10 

11 

1. 

1. 

+EA 

10 

0 

CONM2 

,11, 

11,0,4.986,0. 

,0. 

,0. , ,+EAll 
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+EA11 ,5.79E-02,0.,2. 895E-2 , 0 . ,0. ,2.895E-2 

CONM2, 12 , 11 ,-l, 2. 49E-02, 1.7755 ,-.0572, -.1016, ,+EA12 

+EA12 , .172-05,0. , .55E-05,0. ,0. , .40E-05 

CONM2, 13, 11, -1, 3. 49E-02, 1.7678,0. ,-.1016, ,+EA13 

+EA13 , 3 . 81E-05 , 0 . , .018E-05,0. ,0. , 3.81E-05 

CONH2, 14,11,-1,8.48E-02,1. 7755, . 0572 ,-. 1016, ,+EA14 

+EA14 , . 59E-05 , 0 . ,1.87E-05,0. ,0. ,1.36E-05 

CONM2, 15, 11, -1,3. 49E-02, 1.72 7, 0.0572, -.1016, ,+EA15 

+EA15, . 24E-05 , 0 . , 1 . 07E-05 , 0 . , 0 . , .86E-05 

CONM2, 16, 11,-1,9.47E-02,1.7132, .062, -.1016, ,+EA16 

+E  16,15. 62E-5 , 0 . , .43E-05,0. ,0. ,15.62E-5 

CONM2, 17, 11,-1, 2. 49E-02, 1.7755, .0572, .1016, ,+EA17 

+EA17, .172-05,0., . 55E-05, 0 . , 0. , .40E-05 

CONM2, 18, 11,-1, 3.49E-02, 1.7678, 0. , . 1016, ,+EA18 

+EA18, 3 . 61E-05 , 0 . , . 18E-06, 0 . , 0 . ,3.81E-05 

CONM2 ,19,11,-1,8. 48E-02 ,1.7755, -.0572, .1016, ,+EA19 

+EA19, . 59E-05 , 0  , 1 . 87E-05 , 0 . , 0 . , 1 . 36E-05 

CONM2 ,20, 11 , -1 , 3 . 49E-02 ,1.727, -.0572, .1016, ,+EA20 

+EA20 , . 24E-05 , 0 . ,1.07E-05,0. ,0. , .86E-05 

CONM2 ,21,11,-1,9. 47E-02 , 1.7132, -.062, .1016, ,+EA21 

+EA21, 15. S2E-5,0. , . 43E-05 , 0 . , 0 . , 15 . 62E-5 

CONM2, 22, 11,-1, 51. 36E-2, 1.7678,0. ,-.1016, ,+EA22 

+EA22 , 40 . 99E-5 , 0. , 33 . 13E-5 , 0. , 0. , 40. 99E-5 

CONM2 ,23,11,-1.51. 36E-2 ,1.7678,0.  ,  .U16,  ,+EA23 

+EA23 ,40.99E-5,0. ,33.^3E-5,0. ,0. ,40.99E-5 

MAT1  17 . 42E+102 . 83E+10  .317  2 766 . 3 . 61E-06 

+MA  1 

+MA  1  1.5E+09  1.5E+09  6.8E+07 

PBAR  1  14 . 94E-042 . 71E-081 . 53E-083 . 30E-08  0. 

+PA  1 

+PA  1  0.  0.  0.  0.  0.  0.  0. 

O.+PB  1 

+PB  16.68E-011. 19E+00  0 

SPC  1  1  123456  0. 

SPC  1  2  135  0. 

SPC  1  3  135  0. 

SPC  1  4  135  0. 

SPC  1  5  135  0. 

SPC  1  6  135  0. 

SPC  1  7  135  0. 

SPC  1  8  135  0. 

SPC  1  9  135  0. 

SPC  1  10  135  0. 

SPC  1  11  135  0. 

TLOAD1 ,40,41, ,0,42 
DAREA, 41,11,2, 1000 , 

TABLED1 , 42 ,,,,,,, , +E42 

+E42 , 0 . ,0. ,0.001,0. , .002,1., .003,0. ,+E43 

+E43, .004,0. , ENDT 

TABDMP1 , 99 , CRIT, ,,,,,, +EA99 

+EA99, .1, .065,30. , .023,100., .01, ENDT 

PARAM  AUTOS PC  YES 

EIGR, 97 , MGIV, 0 . , 150. ,3,3, , ,ABC 

+BC , MASS 

TSTEP,98, 10, .001,1,,,, ,+E98 
+E98, , 100, .01,1 
PAR  ul,  POST,  -2 
ENDDATA 
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Appendix  D:  MATLAB  .M  File 


This  MATLAB  .M  file  was  written  to  perform  the  control  optimization 
calculations  in  Chapter  V. 


MATLAB  .M  File 


%This  program  uses  the  MATLAB  subroutine  LQR  to  solve 
%for  the  gain  K  where  u  =  -Kx.  In  addition,  the  Riccati 
%matrix  is  given.  With  the  K  &  S  matrices,  the  cost  functions 
%J  and  C  are  found.  This  ".m"  file  supports  three  mode  shapes. 

% 

%The  state-space  equations  ares 
% 

%  x ' =  Ax  +  Bu 

%  y  =  Cx  +  Du 

% 

%The  closed  loop  A  matrix  becomes: 

% 

%  ACL  =  (A  -  B*K) 

% 

%Initial  conditions  can  be  obtained  from  a  finite  element 
%program.  wl,  w2 ,  w3,  si,  s2 ,  s3,  qOl,  q02,  P,  and  k  must  be 
%defined  prior  to  running  this  .m  file.  w#  are  the  natural 
%frequencies  in  rad/sec;  s#  are  the  damping  factors;  qOl  is  the  3x1 
%matrix  of  initial  displacements;  q02  is  the  3x1  matrix  of 
%initial  velocities;  P  is  the  3x3  matrix  of  eigenvectors;  and  k  is 
%the  modal  stiffness  matrix  such  that  PI*K*P=k  for  PI*M*P=m=I. 

% 

diary  on 
% 

%Equation  parameters  are: 
wl 
w2 
w3 
si 
s2 
s3 
P 
k 
% 

%Calculate  the  inverse  of  P  and  k. 

PI=inv(P) 

ki=inv(k) 

% 

%Initial  conditions  given  are. 

qOl 

q02 

%Transformation  of  initial  conditions  to  modal  coordinates. 

n01=PI*q01 

n02=PI*q02 

x0= [ nOl ; n02 ] 

% 

%Define  undamped  A  and  damped  AD  matrix. 
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A= [ 0  00100 
0  0  0  0  1  0 
0  0  0  0  0  1 
-wl^2  00000 
0  -w2^2  0000 
0  0  -w3~2  000] 

% 

% 

AD= [ 0  00100 
0  0  0  0  1  0 
0  0  0  0  0  1 

-wl*2  0  0  -2*wl*sl  0  0 
0  -w2A2  0  0  -2*w2*s2  0 
0  0  -w3*2  0  0  -2*w3xs3 ] 

% 

%Define  the  B  matrix. 
b= [ -1  1  0; 1  1  0;0  0  1] 
bb=P ' *b 

B= [ 0  0  0 ; 0  0  0 ; 0  0  0;bb] 

% 

%Define  the  C  matrix. 

CP=[0  0  0 ; 0  0  0; 0  0  0] 

CV= [ -1  1  0 ; 1  1  0 ; 0  0  1] 

CPP=CP*P 

CVP=CV*P 

C= [ CPP ( 1 , : )  CVP(1, : ) ;CPP(2,  s )  CVP(2, : ) ;CPP(3,  :  )  CVP(3,:)] 
% 

%Define  the  D  matrix. 

D=[0  0  0; 0  0  0 ; 0  0  0] 

% 

%The  gain  matrix  is  found  from  the  Ricatti  equation  by 
%minimizing  the  cost  function 
% 

%  J  -  int (x ' *Q*x  +  u ' *R*u ) dt 

% 

%The  control  theory  weighing  matrices  are  defined  by: 

% 

Q0= [ wl"2  00000 
0  w2A2  0000 
0  0  w3~2  000 

0  0  0  1  0  0 

0  0  0  0  1  0 

0  0  0  0  0  1] 

% 

Q3=0. 3*Q0 
% 

R0=b ' *P*ki*P ' *b 
R3=0.3*R0 
% 

QI0= [ 150  0  0  000 

0  6750  0  000 

0  0  30375  000 

0  0  0  1  0  0 

0  0  0  0  1  0 

0  0  0  0  0  1] 

% 

QI3=0. 3*QI0 
% 

RI0= ( 0 . 003  0  0 ; 0  0.003  0;0  0  0.0001] 

RI3=0. 3*RI0 
% 

% 
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%Use  MATLAB  subroutine  "lqr”  to  solve  the  Riccati  equation 
%for  the  Riccati  matrix  S  and  the  gain  matrix  K. 

% 

[K00 , S00 ] =lqr { A, B  ,Q0, RO ) 

[K30,S30l=lqr(A,B,Q3,R0) 

[K03 , S03 ] =lqr (A, B,QO , R3 ) 

[K33,S33 ]=lqr (A,B,Q3,R3) 

[ KDOO , SDOO ] = lqr ( AD , B , QO , RO ) 

[KD30, SD30 ] =lqr (AD, B,Q3,R0) 

[ KD03 , SD03 ] =lqr { AD , B,QO, R3 ) 

[KD33 , SD33 ] =lqr (AD, B,Q3 , R3 ) 

[ KIOO, SIOO] =lqr (A, B,QIO, RIO) 

[ KI30 , SI30 ] =lqr ( A, B,QI3 , RIO) 

[ KI03 , SI03 ] =lqr ( A, B,QIO, RI3 ) 

[KI33 , SI33 ] =lqr (A,B,QI3,RI3) 

[ KIDOO, SIDOO ] =lqr (AD, B,QI0 , RIO ) 

[KID30, SID30]=lqr (AD, B,QI3,RI0) 

[ KID03 , SID03 ] =lqr ( AD , B,QIO, RI3 ) 

[KID33,SID33] =lqr ( AD , B, QI3 , RI3 ) 

% 

%Initial  conditions  are  used  to  determine  the  minimum  J 
%values  of  the  control  cost  function. 

% 

JOO=xO ' *SOO*xO 
J30=x0 ' *S30*x0 
J03=x0 ' *S03*x0 
J33=x0 ' *S33*xO 
JDOO=xO ' *SDOO*xO 
JD30=x0 ' *SD30*x0 
JD03=x0 ' *SD03*x0 
JD33=xO ' *SD33*xO 
JIOO=xO ' *SIOO*xO 
JI30=x0' *SI30*x0 
Jl03=x0 ' *SI03*x0 
JI33=xO ' *Sl33*xO 
JIDOO=xO ' *SIDOO*xO 
JID30=x0 ' *SID30*x0 
JID03=x0 ' *SID03*x0 
JID33=xO ' *SID33*xO 
% 

%Initial  control  inputs  are  determined  from  the  control  law. 
% 

FOO=-l*KOO*xO 

F30=-l*K30*x0 

F03=-l*K03*x0 

F33=-l*K33*xO 

FDOO=-l*KDOO*xO 

FD30=-l*KD30*x0 

FD03=-l*KD03*x0 

FD33=-l*KD33*xO 

FIOO=-l*KIOO*xO 

FI30=-l*Kl30*x0 

FI03=-l*KI03*x0 

FI33=-l*KI33*xO 

FIDOO=-l*KIDOO*xO 

FID30=-l*KID30*x0 

FID03=-l*KID03*x0 

FID33=-l*KID33*xO 

% 
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%Expected  values  of  the  control  cost  function  C  =  u'*R*u 
%are  determined. 

% 

COO=FOO ' *RO*FOO 
C30=F30 ' *R0*F30 
C03=F03 ' *R3*F03 
C33=F33’*R3*F33 
CDOO=FDOO ' *RO*FDOO 
CD30=FD30  *  *R0*FD30 
CD03=FD03 ' *R3*FD03 
CD33=FD33 ’ *R3*FD33 
CIOO=FIOO ' *RI0*FI00 
CI30=FI30 '  *RI0*FT30 
CI03=FI03 *  *RI3*FI03 
CI33=FI33 ' *RI3*FI33 
CIDOO=FIDOO ' *RIO*FIDOO 
CID30=FID30 ' *RI0*FID30 
CID03=FID03 ’ *RI3*FID03 
CID33=FID33 ’ *RI3*FID33 
% 

%Open  loop  frequencies  and  damping  factors  are  calculated. 

% 

[  WA, ZA] =damp( A) 

[ WAD , ZAD ] =damp ( AD ) 

% 

%Closed  loop  frequencies  and  damping  factors  are  calculated. 
% 

[WOO, ZOO ] =damp(A-B*KOO ) 

(W30, Z30]=damp(A-B*K30) 

[ W03 , Z03 ] =damp(A-B*K03 ) 

[ W33 , Z33 ] =damp( A-B*K33 ) 

[ WDOO, ZDOO ] =damp( AD-B*KD00 ) 

[ WD30 , ZD30 ] =damp ( AD-B*KD30 ) 

[ WD03 , ZD03 ] =damp( AD-B*KD03 ) 

[ WD33 , ZD33 ] =damp( AD-B*KD33 ) 

[WIOO, ZIOO]=damp(A-B*KlOO) 

[WI30, ZI30]=damp(A-B*Kl30) 

[WI03 , ZI03 ]=damp(A-B*KI03 ) 

[ WI33 , ZI33 ] =damp( A-B*KI33 ) 

[ WIDOO , ZIDOO ] =damp ( AD-B*KIDOO ) 

[WID30, ZlD30]=damp(AD-B*KID30) 

[ WID03 , ZID03 ] =damp( AD-B*KID03 ) 

[ WID33 , ZID33 ] =damp(AD-B*KID33 ) 

% 

diary  off 
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Appendix  E 


and  R  Weights 


Table  E.l  Q  and  R  Nomenclature 


150  0  0  0  0  0 

0  6750  0  000 

0  0  30375  0  0  0 

0  0  0  1  0  0 

0  0  0  0  1  0 

0  0  0  0  0  1 

0.003  0  0 

0  0.003  0 

0  0  0.0001 


(E.l) 


(E.2) 


E-l 


wl 


12.1440 

w2  = 

82.1190 

w3  = 

174.2800 

si  = 

0.0640 

s2  = 

0.0084 

s3  = 

0.0236 

P  = 

0.3759  0  0.0823 

0  -3.6022  0 

0.0775  0  -0.6853 

k  = 

1.0e+004  * 


0.0147 

0 

0 

0 

0.6744 

0 

0 

0 

3.0374 

qOl  = 

1.0e-003  * 

0.1480 

0 

-0.0268 


F-l 


0.1479 

0 

-0.0263 


xO  = 

0.0004 

0 

0.0001 

0.3758 

0 

0.0808 


1.0e+004  * 


0 

0 

0 

0.0001 

0 

0 

0 

0 

0 

0 

0.0001 

0 

0 

0 

0 

0 

0 

0.0001 

-0.0147 

0 

0 

0 

0 

0 

0 

-0.6744 

0 

0 

0 

0 

0 

0 

-3.0374 

0 

0 

0 

AD  = 

1.0e+004  * 


0 

0 

0 

0.0001 

0 

0 

0 

0 

0 

0 

0.0001 

0 

0 

0 

0 

0 

0 

0.0001 

-0.0147 

0 

0 

-0.0002 

0 

0 

0 

-0.6744 

0 

0 

-0.0001 

0 

0 

0 

-3.0374 

0 

0 

-0.0008 

-110 
110 
0  0  1 


B  = 


0 

0 

0 

-0.3759 

-3.6022 

-0.0823 


0 

0 

0 

0.3759 

-3.6022 

0.0823 


0 

0 

0 

0.0775 

0 

-0.6853 


cv  = 


-1 

1 

0 

1 

1 

0 

0 

0 

1 

c  = 

0 

0 

0 

0 

0 

0 

0 

0 

0 

D  = 

0 

0 

0 

0 

0 

0 

0 

0 

0 

Q3  = 

1.0e+003 

★ 

0.0442 

0 

0 

0 

2.0231 

0 

0 

0 

9.1121 

0 

0 

0 

0 

0 

0 

0 

0 

0 

R0  = 

0.0029 

0.0010 

-0.0002 

0.0010 

0.0029 

0.0002 

-0.0002 

0.0002 

0.0001 

KD30  = 

1.0e+003 

★ 

-0.0268 

-0.1312 

-0.6254 

0.0268 

-0.1312 

0.6254 

0.0064 

0.0000 

-6.0627 

SD30  = 

9.1418 

-0.0000 

-0.0000 

0.0000 

69.9676 

-0.0000 

-0.0000 

-0.0000 

143.4420 

0.1402 

0.0000 

0.0000 

0.0000 

0.1402 

-0.0000 

0.0000 

-0.0000 

0.1402 

-0.3759 

-3.6022 

-0.0823 

0.3759 

-3.6022 

0.0823 

0.0775 

0 

-0.6853 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0.0003 

0 

0 

0 

0.0003 

0 

0 

0 

0.0003 

-0.0102 

-0.0085 

-0.0183 

0.0102 

-0.0085 

0.0183 

0.0024 

0.0000 

-0.1777 

0.1402 

-0.0000 

-0.0000 

-0.0000 

0.1402 

0.0000 

-0.0000 

-0.0000 

0.1402 

0.0531 

0.0000 

0.0000 

-0.0000 

0.0091 

0.0000 

-0.0000 

-0.0000 

0.0041 

F-3 


JD30 


0.0076 


FD30  = 


5.3613 

-5.3613 

13.9415 


CD30  = 


0.0626 


WD30  = 

12.9672 

12.9672 

186.0946 

186.0946 

87.6858 

87.6858 


ZD30  = 


0.3617 

0.3617 

0.3574 

0.3574 

0.3568 

0.3568 
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